一階常微分方程式をホイン法と求積法を用いたグラフの描き方を教えてください。

[AMG]

この問題は例のどの部分をどう変えれば良いのかを教えていただけたら幸いです。C言語に関してはあまり詳しくないのでかみ砕いて説明していただけると大変ありがたいです。
この掲示板を初めて使用させていただいているので、質問の仕方が間違っているなど、失礼があればぜひ教えて下さい。
よろしくお願いいたします。

次の一階常微分方程式（ロジスティック方程式）
dx/dt=x-x^2 　　…(1)式
をt= [0,20]で初期条件 x(0) = 0.1のもとにホイン法を適用して数値的に解いてグラフを描き、求積法のよって求めた式を用いて描いたx(t)のグラフと比較せよ．ただし差分間隔 ∆ を0.1 とする．x(20) の求積法による解は 0.999999981，ホイン法による数値解も 0.999999981 となり，正確に計算できていることが分かる．ロジスティク方程式はいろいろな物理現象を表しているが，この方程式の差分近似式とカオスの関係について調べよ。

例として
dx/dt=x
について各方法の精度を比較してみよう．この一階常微分方程式を前進オイラー法，ホイン法および４次のルンゲ・クッタ法を用いてt=[0,1]の範囲で数値計算し，それぞれの方法で求めたx(1) の値に含まれる誤差はどのようになるか検討する．ただし，初期条件をx(0) = 0.1 ，差分間隔 ∆ を 0.01 とする． まず解析解を求めると，

∫(dx/x)=∫dt
lnx =t + C
x=e^(t+C) = C′*e^t
初期条件より
x(0) = 0.1 = C′e^0 =C′ ⟺ C′= 0.1
∴ x(t) = 0.1e^t x(1) =e/10

オイラー法を用いると，
f(xi,ti) =xi
xi+1= xi+f(xi,ti)∆t
で，プログラムの例を以下に示す．
dt = 0.01;
n = (int) (1.0/0.01);
x[0] = 0.1; t[0] = 0.0;
for (i = 1; i <= n; i++)
{
x = x[i-1] + func[x[i-1]] * dt;
t = (double) i * dt;
}
double func(double x)
{
return x;
}

ホイン法を用いると，
f(xi,ti) =xi
k1= f(xi,ti)∆t
k2= f(xi+k1,ti+∆t)∆t
xi+1=xi+1/2(k1+k2)
で，プログラムの例を以下に示す．
dt = 0.01;
n = (int) (1.0/0.01);
x[0] = 0.1; t[0] = 0.0;
for (i = 1; i < = n; i++)
{
k1 =func(x[i-1]) * dt;
k2 = func(x[i-1] + k1) * dt;
x = x[i-1] + 0.5 * (k1 + k2);
t = (double) i * dt;
}
double func(double x)
{
return x;
}


（参考）求積法による解
dx/dt=ax-bx^2 (a≠0)
⇔dx/(ax-bx^2)=dt
⇔dx/(x(1-bx/a))=adt
⇔dx/x+(b/a*dx)/(1-b/a*x)=adt
1-b/a*x=yとおくと、-b/a*dx=d
　　∴dx/x-dy/y=adt
⇔∫(dx/x)-∫(dy/y)=a∫dt
⇔Inx-Iny=at+C
⇔In(x/y)=at+C
⇔x=yexp(at+C)
⇔x=(1-bx/a)exp(at+C)=C(a-bx/a)exp(at)
⇔x=(Caexp(at))/(a+Cbexp(at))=Ca/(aexp(-at)+Cb)
x(0)=x0とすると、
　　　　　x0=x(0)=Ca/(aexp(-a0)+Cb)=Ca/a+Cb
∴C=ax0/(a-bx0)

したがって、
x(t)=a(ax0/(a-bx0)) / (aexp(-at)+b(ax0/(a-bx0)))
=(a^2)x0/(a(a-bx0)exp(-at)+abx0)
(1)式では、a=1、b=1、また初期値としてx0=0.1なので、
x(t)=0.1/(0.9exp(-20)+0.1)=0.999999981
である。

(参考)一階常微分方程式の解法
次の一階常微分方程式
　　　　dx(t)/dt=f(x,t) …(2)式
を用いて初期条件x(0) =x0からx(t)の値を順に数値的に求めることを考える．その際区間[0,T]をN等分して∆t =T/Nとし，ti=i∆t(i= 0,1,2,⋯,N) とおき,t=[0,T]内の連続した時刻におけるx(t) を考えるのではなく，とびとびのx(ti)のみを考える．

オイラー法
dx/dtの前進差分近似式は，導関数の定義より，
dx(t)/dt = lim ∆t→0 (x(t+Δt)-x(t))/Δt≈(x(t+Δt)-x(t))/Δt
と表せる．したがって(2)式は,
dx(t)/dt=f(x,t)≈(x(t+Δt)-x(t))/Δt
⇔x(t+Δt)≈x(t)+f(x,t)Δt
となる．また離散点に関する式で表わすと，
xi+1=xi+f(xi,ti)Δt ti=iΔt i=0,1,2,…,N-1
となるので，x0からx1を,x1からx2を求め，計算を繰り返してxN=x(T) まで求めることができる．

ホイン法
２次のルンゲ・クッタ法の代表的なものにホイン法がある．一般に２次のルンゲ・クッタ法は以下の式で表わされる．

k1=f(x,t)Δt
k2=f(x+pk1,t+qΔt)Δt
x(t+Δt)=x(t)+ak1+bk2 …3式合わせて(3)式
ここで，
f(x+pk1,t+qΔt)=f(x,t)+∂f/∂x *pk1+∂f/∂t *qΔt+o(Δt^2)
と展開したものを(3)式のk2に代入すると，
x(t+Δt)=x(t)+af(x,t)Δt+b{f(x,t)+∂f/∂x *pfΔt+∂f/∂t *qΔt+0(Δt^2)}Δt
　　　　⇔x(t+Δt)=x(t)+(a+b)f(x,t)Δt+bpf *∂f/∂xΔt^2+bq *∂f/∂t *Δt^2+o(Δt^3)
　　　　　　　…(4)式
(4)式から，２次のルンゲ・クッタ法がo(∆t^3)の誤差を持つ，すなわち２次の精度を持つことは明らかである. x(t+∆t)のテイラー展開は先に示したように，
x(t+∆t) = x(t)+x′(t)∆t+1/2*(x^2)x′′(t)∆t^2 +o(∆t^3)
と表せる．この式に，
x′(t) = f(x,t)
x′′(t) =∂f/∂x*∂x/∂t+∂f/∂t*∂t/∂t
=∂f/∂x*f+∂f/∂t
を代入すると，
x(t+∆t) =x(t)+f(x,t)∆t+1/2(∂f/∂x *f+∂f/∂t)Δt^2+o(∆t^3)
=x(t)+f(x,t)∆t+1/2 *f*∂f/∂xΔt^2+1/2*∂f/∂t∆t^2+o(Δt^3)
…(5)式
(4)式と(5)式の比較から，
a+b= 1, bp=1/2, bq=1/2
となり，1つの解は
a=b=1/2， p=q= 1
である．これを1.10.5式に代入すると，
k1=f(xi,ti)Δt
k2= f(xi+k1,ti+Δt)Δt
x(t+Δt)=x(t)+1/2*(k1+k2)
となる．これがホイン法です。

[みけCAT]

課題の丸投げは禁止です。
フォーラムルールをお読みください。
AMGさんはどこまでわかっていて、何がわからない(聞きたい)のかを具体的に質問していただけると幸いです。

「この問題は例のどの部分をどう変えれば良いのか」という質問が提示されているので、よく見たらそれほど悪い投稿ではない気がした。

[みけCAT]

コード:

dt = 0.01; /* 差分間隔 */
n = (int) (1.0/0.01); /* 区間の終わり÷終端間隔 */
/* ※終端間隔はリテラルではなくdtを使うべきである。
     また、単純に実数の割り算をして整数に変換すると、誤差による誤動作が心配である。 */
x[0] = 0.1; t[0] = 0.0; /* 関数の値(x)とその引数(t)の初期値 */
 for (i = 1; i < = n; i++)
    {
       k1 =func(x[i-1]) * dt;
       k2 = func(x[i-1] + k1) * dt;
       x[i] = x[i-1] + 0.5 * (k1 + k2);
       t[i] = (double) i * dt;
    }
 double func(double x)
 {
 return x; /* 対象の関数のdx/dt */
 }

このへんが解析対象のパラメータのようなので、問題に合わせて変えるといいでしょう。
また、ループの条件はi < = nではなく、i <= n (<演算子と=演算子を並べるのではなく、<=演算子を用いる)とした方がいいでしょう。
さらに、配列xおよびtの要素数を十分確保しないといけません。