I'm trying to solve this numerical differential equations, can someone help?
clc;
clear;
syms A1(z) A2(z)
lamda1 = 1560*(10^-9);
c=3*(10^8);
d_eff=27*(10^-12);
omga1=(2*pi*c)/(lamda1);
omga2=omga1*2;
n=2.2;
k1=(n*omga1)/c;
k2=(n*omga2)/c;
ode1 = diff(A1) == (2*i*(omga1^2)*d_eff*A2*conj(A1)*exp(-i*(2*k1-k2)*z))/(k1*(c^2));
ode2 = diff(A2) == (i*(omga2^2)*d_eff.*(A1.^2).*exp(i*(2*k1-k2)*z))/(k2*(c^2));
odes = [ode1; ode2];
cond1 = A1(0) == 1;
cond2 = A2(0) == 0;
conds = [cond1; cond2];
M = matlabFunction(odes)
sol = ode45(M(1),[0 20],[2 0]);
in this question both ODE are coupled, hence there's only 1 ODE to solve:
1.- use 1st equation to write
and feed this expression into 2nd equation.
2.- regroup
now the ODE to solve is
3.- It can obviously be done in MATLAB, but for this particular ODE in my opinion the Wolfram online ODE Symbolic Solver does a better job.
Input the obtained ODE of previous point into the ODE solver available in this link
https://www.wolframalpha.com/input?i=y%27%27+%2B+y+%3D+0
and solve
4.- The general (symbolic) solutions for
A2
areNote that I used
k1
instead ofn1
andk2
instead ofn2
just in the Wolfram ODE Solver.Rewording ; the
k1
k2
expressions of the general solutions are not the wave numbersk1
k2
of the equations in the question. Just replace accordingly.5.- Now get
A1
using expression in point 2.6.- I have spotted 2 possible errors in the MATLAB code posted in the question that shouldn't be ignored by the question originator:
In the far right side of the posted MATLAB code, the exponential expressions
6.1.- both show
exp(-i(2*k1-k2)z))/(k1(c^2))
there's this
(k1*(c^2))
dividing the exponent wheras in the question none of the exponentials show such denominator their respective exponents.6.2.- the
dk
ordelta k
expression in the exponentials of the question are obviouslyk2-k1
ork1-k2
, here there may be room for a sign ambiguity, that may shoot a wave solution onto the opposite direction, yet the point here is whereshould probably be
or just
6.3.- and yes, in MATLAB
(-1)^.5
can be either expressed with1j
or1i
but as written in the MATLAB code made available in the question, since only a chunk of code has been made available, it's fair to assume that no suchi=1j
has been done.