Error in ode45 solution

3 次查看（过去 30 天）
PS
PS 2024-5-30
编辑： David Goodmanson 2024-6-2
Ran in:
Ode45 results dont maintain power conservation for some cases. How to ensure that power conservation is maintained?
While using ode45 to solve coupled differential equation, when one of the input parameters is real,the solution to ode45 maintains power conservation. When i change that input parameter from real to complex then ode45 solution is not maintaing power conservation. I have used relatve error tolerance = 1e-6 and absolute error tolerance = 1e-9. Please let me know what can be done to get the correct results.
For eg, In the below code, tv* parameter was initially real, which maintained power conservation, however when i change tv* to complex, the solution doesnt maintain power conservation.
options=odeset('RelTol',1e-6,'AbsTol',1e-9);
f=@(z,xx_val) -1i*[tv12*xx_val(2)*exp(1i*del_beta12*z)+tv13*xx_val(3)*exp(1i*del_beta13*z)+tv14*xx_val(4)*exp(1i*del_beta14*z)+tv15*xx_val(5)*exp(1i*del_beta15*z)+tv16*xx_val(6)*exp(1i*del_beta16*z);...
tv21*xx_val(1)*exp(1i*del_beta21*z)+tv23*xx_val(3)*exp(1i*del_beta23*z)+tv24*xx_val(4)*exp(1i*del_beta24*z)+tv25*xx_val(5)*exp(1i*del_beta25*z)+tv26*xx_val(6)*exp(1i*del_beta26*z);...
tv31*xx_val(1)*exp(1i*del_beta31*z)+tv32*xx_val(2)*exp(1i*del_beta32*z)+tv34*xx_val(4)*exp(1i*del_beta34*z)+tv35*xx_val(5)*exp(1i*del_beta35*z)+tv36*xx_val(6)*exp(1i*del_beta36*z);...
tv41*xx_val(1)*exp(1i*del_beta41*z)+tv42*xx_val(2)*exp(1i*del_beta42*z)+tv43*xx_val(3)*exp(1i*del_beta43*z)+tv45*xx_val(5)*exp(1i*del_beta45*z)+tv46*xx_val(6)*exp(1i*del_beta46*z);...
tv51*xx_val(1)*exp(1i*del_beta51*z)+tv52*xx_val(2)*exp(1i*del_beta52*z)+tv53*xx_val(3)*exp(1i*del_beta53*z)+tv54*xx_val(4)*exp(1i*del_beta54*z)+tv56*xx_val(6)*exp(1i*del_beta56*z);...
tv61*xx_val(1)*exp(1i*del_beta61*z)+tv62*xx_val(2)*exp(1i*del_beta62*z)+tv63*xx_val(3)*exp(1i*del_beta63*z)+tv64*xx_val(4)*exp(1i*del_beta64*z)+tv65*xx_val(5)*exp(1i*del_beta65*z)];
[zz,xa_var]=ode45(f,[0 1],init_condit(nn,:),options); %[LP01 LP11a] represent the initial condition
Unrecognized function or variable 'init_condit'.
  2 个评论
Torsten
Torsten 2024-5-30
编辑：Torsten 2024-5-30
Your code should be execuable and reproduce the problem you are describing. Further - for us uneducated forum people - you will need to explain how energy conservation can be checked from the results you get.
David Goodmanson
David Goodmanson 2024-6-1
编辑：David Goodmanson 2024-6-1
Hi PS, is it true that (before you start changing things around) all of the del_beta variables are real just like all of the tv variables are, and that for example tv23 = tv32 and del_beta23 = -del_beta32, and the same for the others?

请先登录，再进行评论。

请先登录，再回答此问题。

回答（1 个）

David Goodmanson
David Goodmanson 2024-6-2
编辑：David Goodmanson 2024-6-2
Hello PS,
Assuming the following (before variables get changed around)
tv_jk = tv_kj, for all j,k all real
del_beta_jk = -del_beta_kj for all j,k all real
and that for the 6x1 solution vector xx_val = y, what you refer to as power is
power = <y|y> = norm(y)^2 = y'*y = sum(y.*conj(y))
Then:
The equation is
i*dydt = H*y
and with the assumptions on variables above, H is Hermitian, H = H' and consequently
<y|y> is constant in time.
I am not sure what you mean by 'one of the input parameters' not being real, but if parameters are changed in H to make it not Hermitian, then in general <y|y> will not be conserved. (There may be some odd exception where it's still conserved, possibly). But if you change some of the tv to complex and keep with the condition
tv_jk = conj(tv_kj) % for all j,k
then H remains Hermitian and <y|y> should be constant in time.

类别

在 Help Center 和 File Exchange 中查找有关 Ordinary Differential Equations 的更多信息

标签

产品


版本

R2019a

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by