[Date Index]
[Thread Index]
[Author Index]
Re: DSolve for a real function
*To*: mathgroup at smc.vnet.net
*Subject*: [mg127942] Re: DSolve for a real function
*From*: Andreas Talmon l'Armée at smc.vnet.net
*Date*: Tue, 4 Sep 2012 05:46:06 -0400 (EDT)
*Delivered-to*: l-mathgroup@mail-archive0.wolfram.com
*Delivered-to*: l-mathgroup@wolfram.com
*Delivered-to*: mathgroup-newout@smc.vnet.net
*Delivered-to*: mathgroup-newsend@smc.vnet.net
*References*: <k1k8np$90f$1@smc.vnet.net> <20120901062746.B3D7F687B@smc.vnet.net> <504462A6.5070404@fsm.tu-darmstadt.de> <CAEtRDSfhfr=mca7_3XsfGrF_5MbAuvAF5f8oVo7y+M9=_C68Ow@mail.gmail.com> <5044BD10.5070102@fsm.tu-darmstadt.de> <CAEtRDSdHoBwvVi3ewxnnnQD_rrT_+PFDKi1sHnWct3L9eJQt0g@mail.gmail.com>
Hi,
I think this is not a numerical problem. Because if I look at the
outcome of the symbolic computations for y3[x] (Simplify) and y4[x]
(FullSimplify) there is one solution without imaginary part (Simplify)
and one with imaginary part (FullSimplify) .
http://dl.dropbox.com/u/4920002/DGL_4th_Order_with_own_solution.nb
I am feeling a bit uncomfortable with these outcomes because I think the
solution y3[x] (Simplify) is the better one even if the real part of
both solutions is the same. But I think I just have to accept it. Is
there a way to suppress the calculation with complex numbers.
Thanks for your time and your help
Andreas Talmon
On 09/03/2012 06:54 PM, Bob Hanlon wrote:
> In numerical computations, the specific computations and their order
> of computation determine the accumulation of the numerical noise. You
> generally reduce the accumulation of numerical noise by increasing the
> working precision.
>
> For example in your t2 rather than using machine precision specify
> extended precision (say 55 places)
>
> t2 = Table[
> y4[x] /. {c -> 1, ic0 -> 2, ic1 -> 3, a -> 4}, {x, 0, 100, 10}] //N[#, 55] &
>
> Similarly, when plotting y4
>
> pRe2 = Plot[
> Chop[y4[x]] /. {c -> 1, ic0 -> 2, ic1 -> 3, a -> 4}, {x, 0, 100},
> WorkingPrecision -> 55, PlotStyle -> Red];
>
> pIm2 = Plot[Im[y4[x] /. {c -> 1, ic0 -> 2, ic1 -> 3, a -> 4}],
> {x, 0, 100}, WorkingPrecision -> 55,
> Frame -> True, Axes -> False]
>
>
> Bob Hanlon
>
>
> On Mon, Sep 3, 2012 at 10:22 AM, Andreas Talmon l'Armée
> <talmon at fsm.tu-darmstadt.de> wrote:
>> You are completely right but I think there should be no real part in the
>> solution.
>>
>> I did the solution of the differential equation on paper and got the term
>> for y1[x].
>>
>> Then I determined the constants C[5]...C[8] with the boundary conditions
>> with Solve[] in Mathematica. There I got an unusual behaviour.
>>
>> When using Simplify on the sol1 could generate a soltuion with only a real
>> part.
>> When using FullSimplify I get the imaginary part which is very small for
>> small x but for greater x also the imaginary part starts to oscillate with a
>> great amplitude.
>>
>> I do not understand how a solution can be Real or Complex depending on the
>> Simplifications scheme I use.
>>
>> Shouldn't there be one solution to a problem with specified boundary
>> conditions.?
>>
>> The notebook where I did all this is on:
>> http://dl.dropbox.com/u/4920002/DGL_4th_Order_with_own_solution.nb
>>
>>
>> Thanks for your great help,
>>
>> Andreas.
>>
>>
>>
>>
>> On 09/03/2012 03:22 PM, Bob Hanlon wrote:
>>> Clear[sol];
>>>
>>> $Assumptions = Element[{a, c, ic0, ic1}, Reals];
>>>
>>> sol[a_, x_] =
>>> y[x] /. DSolve[{y''''[x] + a y[x] == 0, y''[-c] == ic0, y''[c] == ic0,
>>> y'''[-c] == ic1, y'''[c] == -ic1}, y[x], x][[1]] // FullSimplify
>>>
>>> ((1 + I)*a^(1/4)*
>>> (I*E^(I*Sqrt[2]*a^(1/4)*c) -
>>> I*E^(Sqrt[2]*a^(1/4)*c) +
>>> E^(Sqrt[2]*a^(1/4)*
>>> ((1 + I)*c + I*x)) -
>>> E^(I*Sqrt[2]*a^(1/4)*x) -
>>> E^(Sqrt[2]*a^(1/4)*x) +
>>> E^(Sqrt[2]*a^(1/4)*
>>> ((1 + I)*c + x)) +
>>> I*E^(Sqrt[2]*a^(1/4)*
>>> (I*c + (1 + I)*x)) -
>>> I*E^(Sqrt[2]*a^(1/4)*
>>> (c + (1 + I)*x)))*ic0 +
>>> Sqrt[2]*(E^(I*Sqrt[2]*a^(1/4)*
>>> c) + E^(Sqrt[2]*a^(1/4)*
>>> c) + E^(Sqrt[2]*a^(1/4)*
>>> ((1 + I)*c + I*x)) +
>>> E^(I*Sqrt[2]*a^(1/4)*x) +
>>> E^(Sqrt[2]*a^(1/4)*x) +
>>> E^(Sqrt[2]*a^(1/4)*
>>> ((1 + I)*c + x)) +
>>> E^(Sqrt[2]*a^(1/4)*
>>> (I*c + (1 + I)*x)) +
>>> E^(Sqrt[2]*a^(1/4)*
>>> (c + (1 + I)*x)))*ic1)/
>>> E^((-1)^(1/4)*a^(1/4)*(c + x))/
>>> (4*a^(3/4)*(Sin[Sqrt[2]*a^(1/4)*
>>> c] + Sinh[Sqrt[2]*a^(1/4)*
>>> c]))
>>>
>>> Table[sol[a, x] /. {
>>> c -> RandomReal[{0, 10}, WorkingPrecision -> 20],
>>> ic0 -> RandomReal[{0, 10}, WorkingPrecision -> 20],
>>> ic1 -> RandomReal[{0, 10}, WorkingPrecision -> 20]},
>>> {a, 1, 5}, {x, 0, 4}]
>>>
>>> {{0.01540650723811880994 + 0.*10^-21 I, -0.0003703423761496348 +
>>> 0.*10^-20 I, -0.1278269715054133646 + 0.*10^-20 I,
>>> 0.0392723836429585004 + 0.*10^-20 I,
>>> 13.10841701158527505 + 0.*10^-18 I}, {-0.00286328430472954469 +
>>> 0.*10^-21 I, -1.985775360248369253 + 0.*10^-19 I,
>>> -0.311655087881410281 +
>>> 0.*10^-19 I, -0.718055887186593765 + 0.*10^-19 I,
>>> 21.42291506291193811 + 0.*10^-18 I}, {0.419394450195258310 + 0.*10^-19
>>> I,
>>> 0.01666534839706064001 + 0.*10^-21 I,
>>> 0.528881084742103412 + 0.*10^-19 I, -0.060937975960727530 + 0.*10^-19
>>> I,
>>> 0.0086995512380544778 + 0.*10^-20 I}, {0.01847749159520439918 +
>>> 0.*10^-21 I,
>>> 0.00517593237021928535 + 0.*10^-21 I,
>>> 0.0134376958759010127 + 0.*10^-20 I, -0.0679397801010012506 +
>>> 0.*10^-20 I, -5.21491473427110550 +
>>> 0.*10^-18 I}, {-0.00168718272676310080 +
>>> 0.*10^-21 I, -0.00137069152420374693 + 0.*10^-21 I,
>>> 3.051457492528748582 + 0.*10^-19 I,
>>> 0.00026536991493989622 + 0.*10^-21 I, -59.38147365750750826 + 0.*10^-18
>>> I}}
>>>
>>> The imaginary parts cancel out; the residual imaginary parts are just
>>> numerical noise. This can be removed with Chop
>>>
>>> % // Chop
>>>
>>> {{0.01540650723811880994, -0.0003703423761496348, -0.1278269715054133646,
>>> 0.0392723836429585004,
>>> 13.10841701158527505}, {-0.00286328430472954469, -1.985775360248369253,
>>> \
>>> -0.311655087881410281, -0.718055887186593765,
>>> 21.42291506291193811}, {0.419394450195258310, 0.01666534839706064001,
>>> 0.528881084742103412, -0.060937975960727530,
>>> 0.0086995512380544778}, {0.01847749159520439918,
>>> 0.00517593237021928535,
>>> 0.0134376958759010127, -0.0679397801010012506, -5.21491473427110550}, \
>>> {-0.00168718272676310080, -0.00137069152420374693, 3.051457492528748582,
>>> 0.00026536991493989622, -59.38147365750750826}}
>>>
>>> ?Chop
>>>
>>> Chop[expr] replaces approximate real numbers in expr that are close to
>>> zero by the exact integer 0. >>
>>>
>>>
>>> Bob Hanlon
>>>
>>>
>>> On Mon, Sep 3, 2012 at 3:56 AM, Andreas Talmon l'Armée
>>> <talmon at fsm.tu-darmstadt.de> wrote:
>>>> Hi,
>>>>
>>>> My initial conditions are the following and I am pretty sure that my
>>>> solution consists only of a real part.
>>>> The Solution has four eigenvalues and they are complex conjugated. With
>>>> all
>>>> variables and all parameters real numbers I must be able to retrieve a
>>>> real
>>>> solution. But being a mathematica newbie, I do not understand how to do
>>>> it
>>>> with mathematca.
>>>>
>>>> y''[-c] == ic0, y''[c] == ic0, y'''[-c] == ic1, y'''[c] == -ic1
>>>>
>>>> My Notebook is also ready for download at dropbox.com:
>>>> http://dl.dropbox.com/u/4920002/DGL_4th_Order.nb
>>>>
>>>> Clear[sol]
>>>>
>>>> $Assumptions = {a \[Element] Reals, ic0 \[Element] Reals,
>>>> ic1 \[Element] Reals, c \[Element] Reals};
>>>>
>>>> sol[a_, x_] =
>>>> y[x] /. DSolve[{y''''[x] + a y[x] == 0, y''[-c] == ic0,
>>>> y''[c] == ic0, y'''[-c] == ic1, y'''[c] == -ic1}, y[x],
>>>> x][[1]] // FullSimplify
>>>>
>>>> Reduce[Element[sol[a, x], Reals], a, Reals]
>>>>
>>>>
>>>>
>>>> Thanks for your help,
>>>>
>>>> Andreas
>>>>
>>>>
>>>>
>>>>
>>>>
>>>>
>>>>
>>>> On 09/01/2012 08:27 AM, Bob Hanlon wrote:
>>>>> What are your initial conditions?
>>>>>
>>>>> Clear[sol]
>>>>>
>>>>> sol[a_, x_] = y[x] /. DSolve[
>>>>> {y''''[x] + a y[x] == 0,
>>>>> y[0] == ic0, y'[0] == ic1,
>>>>> y''[0] == ic2, y'''[0] == ic3},
>>>>> y[x], x][[1]] // FullSimplify
>>>>>
>>>>> (1/(2*a^(3/4)))*
>>>>> (Cosh[(a^(1/4)*x)/Sqrt[2]]*
>>>>> (2*a^(3/4)*ic0*Cos[(a^(1/4)*x)/
>>>>> Sqrt[2]] + Sqrt[2]*
>>>>> (Sqrt[a]*ic1 + ic3)*
>>>>> Sin[(a^(1/4)*x)/Sqrt[2]]) +
>>>>> (Sqrt[2]*(Sqrt[a]*ic1 - ic3)*
>>>>> Cos[(a^(1/4)*x)/Sqrt[2]] +
>>>>> 2*a^(1/4)*ic2*Sin[(a^(1/4)*x)/
>>>>> Sqrt[2]])*
>>>>> Sinh[(a^(1/4)*x)/Sqrt[2]])
>>>>>
>>>>> Reduce[Element[sol[a, x], Reals], a, Reals]
>>>>>
>>>>> a > 0
>>>>>
>>>>>
>>>>> Bob Hanlon
>>>>>
>>>>>
>>>>> On Fri, Aug 31, 2012 at 3:59 AM, <"Andreas Talmon
>>>>> l'Arm=E9e"@smc.vnet.net> wrote:
>>>>>> Hi All
>>>>>>
>>>>>> Is there a way to tell mathematica to solve only for real solutions. My
>>>>>> differential equation is of the kind
>>>>>>
>>>>>> y''''[x]+a y[x]==0
>>>>>>
>>>>>> a= constant coefficient
>>>>>>
>>>>>> I know that I get 4 komplex eigenvalues which are complex conjungated.
>>>>>> But y[x] is a real function.
>>>>>> Solving this equation with DSolve always gets a complex function y[x].
>>>>>>
>>>>>> Any Ideas.
>>>>>>
>>>>>> Thanks, Andreas
>>>>>>
Prev by Date:
**Re: DSolve for a real function**
Next by Date:
**Re: DSolve for a real function**
Previous by thread:
**Re: DSolve for a real function**
Next by thread:
**Re: DSolve for a real function**
| |