Services & Resources / Wolfram Forums
-----
 /
MathGroup Archive
2004
*January
*February
*March
*April
*May
*June
*July
*August
*September
*October
*November
*December
*Archive Index
*Ask about this page
*Print this page
*Give us feedback
*Sign up for the Wolfram Insider

MathGroup Archive 2004

[Date Index] [Thread Index] [Author Index]

Search the Archive

Re: Iterate initial conditions

  • To: mathgroup at smc.vnet.net
  • Subject: [mg48838] Re: Iterate initial conditions
  • From: "Peter Pein" <petsie at arcor.de>
  • Date: Fri, 18 Jun 2004 02:13:08 -0400 (EDT)
  • References: <cap3en$c95$1@smc.vnet.net>
  • Sender: owner-wri-mathgroup at wolfram.com

"Michael Hohendorf" <se0085 at uni-essen.de> schrieb im Newsbeitrag
news:cap3en$c95$1 at smc.vnet.net...
> Hi,
>
> I have a system of differential equations, wich are results in two
> functions f[x], g[x].
> Now I have to iterate the initial conditons
> (Tau_11[0],Tau_22[0],Tau_33[0]) until the functions f[x] and g[x] are
> equal to other functions Sigma_m[x] and Sigma_u[x].
> Can someone help me to solve my problem?
>
> Michael
>
> solution = NDSolve[{
> Tau_11'[x] == (2 a_T[x] v'[x] η[x] - Exp[Epsilon[x] λ[x](Ï"au_11[x]
> + Tau_22[x] + Tau_33[x])/η[x]] Ï"au_11[x] + 2 a_T[x] v'[x] η[x] λ[x]
(1 -
> ξ)) 1/(a_T[x]v[x]λ[x]),
>
> Tau_22'[x] == (2 a_T[x] η[x] (v[x] s'[x]/s[x])- Exp[Epsilon[x]
> λ[x](Ï"au_11[x] + Tau_22[x] + Tau_33[x])/η[x]] Ï"au_22[x] + 2 a_T[x]
η[x]
> λ[x] (v[x] s'[x]/s[x]) (1 - ξ)) 1/(a_T[x]v[x]λ[x]),
>
> Tau_33'[x] == (2 a_T[x] η[x] (v[x] r'[x]/r[x])- Exp[Epsilon[x]
> λ[x](Ï"au_11[x] + Tau_22[x] + Tau_33[x])/η[x]] Ï"au_33[x] + 2 a_T[x]
η[x]
> λ[x] (v[x] r'[x]/r[x])(1 - ξ)) 1/(a_T[x]v[x]λ[x]),
>
> Tau_12'[x] == (a_T[x] λ[x] v'[x] r'[x] Tau_12[x] (1 - ξ) + a_T[x] λ[x]
> Tau_12[x] (v[x] s'[x]/s[x]) (1 - ξ) - Exp[Epsilon[x] λ[x](Ï"au_11[x]
> + Tau_22[x] + Tau_33[x])/η[x]] Ï"au_12[x]) 1/(a_T[x]v[x]λ[x]),
>
> Tau_11[0]==a , Tau_22[0]==b , Tau_33[0]==c , Tau_12[0]==9810},
> {Tau_11[x],Tau_22[x],Tau_33[x],Tau_12[x]},{x,0,0.42}]
>
I suppose the following should read
Tau_11[x_]=Tau_11[x]/.solution[[1]]
Tau_22[x_]=Tau_22[x]/.solution[[1]]
Tau_33[x_]=Tau_33[x]/.solution[[1]]
Tau_12[x_]=Tau_12[x]/.solution[[1]]
> Tau_11[x_]=Tau_11[x]/.solution[[1]]
> Tau_22[x_]=Tau_11[x]/.solution[[1]]
> Tau_33[x_]=Tau_11[x]/.solution[[1]]
> Tau_12[x_]=Tau_11[x]/.solution[[1]]
>
> f[x_]=Tau_11[x]-Tau_22[x]
> g[x_]=Tau_22[x]-Tau_33[x]
>
> Iterate initial conditions until
>
> f[x]==Sigma_m[x]
> g[x]==Sigma_u[x]
>
Hi Michael,

 as Alois already emphasized, it's sufficient to calculate Sigma_m[0] and
Sigma_u[0] and choose the inital conditions to satisfy
S_m[0]==T_11[0]-T_22[0] and S_m[0]==T_22[0]-T_33[0].
But if the given solutions Sigma_x are not exact combinations of the
obtained solutions Tau_ij, you have to iterate the init.conds until the
error becomes minimal. I'll demonstrate it for an easy example, because your
posting contains a lot of unreadable (resp. meaningless (at least to me))
characters.
If we have got an diff-eq. which "looks almost as" s''+4s=0, we can try to
find a solution which approximates a known solution to the simple DEq:

In[1]:= err[target_, (s0_)?NumericQ, (s1_)?NumericQ] :=
  Module[{sf, s, t},
   sglobal = sf = s /. First[NDSolve[{Derivative[2][s][t] -
Derivative[1][s][t]/
             (100 + t) + 4.01*s[t] == Exp[-(1 + t^2)]*Sin[5*t], s[0] == s0,
          Derivative[1][s][0] == s1}, s, {t, -Pi, Pi}]];
    NIntegrate[(sf[t] - target[t])^2, {t, -Pi, 0, Pi}]]

this is one known solution for s''+4s==0:
In[2]:= test = Cos[2*(#1 - 3)] + Sin[2*(#1 + 1)] & ;

In[3]:= FindMinimum[err[test, s0, s1], {s0, 1, 2}, {s1, -2, -1}]
Out[3]= {0.0011098774260274377, {s0 -> 1.867576215614434,
s1 -> -1.4803450065191435}}

this would have been s0 and s1 for the simple DEq.
In[4]:= N[{test[0], Derivative[1][test][0]}]
Out[4]= {1.8694677134760478, -1.3911246694921364}

In[5]:= Plot[sglobal[t] - test[t], {t, -Pi, Pi}];

You can play around with the NIntegrate-statement, to obtain weighted errors
etc.

Have fun ;-)
-- 
Peter Pein, Berlin
to write to me, start the subject with [




  • Prev by Date: Re: Symbolic use of numerical function FindRoot via ?NumericQ
  • Next by Date: Mathematica can't do this double integral
  • Previous by thread: Re: Iterate initial conditions
  • Next by thread: Suggestions for more compact code