How to solve the continuity equation for dark and light scenario by symbolic function with boundary condition
4 views (last 30 days)
Show older comments
Jiaxing Liu
on 15 Jul 2021
Commented: Jiaxing Liu
on 16 Jul 2021
Now I met a problem on my research. I met with two ODEs with identical form in left side but one is homogenous and the other is nonhomogeneous.
Here is my code:
syms V_T mu_p F r_m % constant of ODE1
syms N_v phi_3 phi_4 d_j % constants for boundary condition
syms p(x) G(x) % p means hole density G(x) means generation rate
eqn = -V_T*mu_p*diff(p,x,2) + F*mu_p*diff(p,x) + r_m*p == 0; %MARKED
P_dark = dsolve(eqn) %% no ; so the form is readable
eqn = -V_T*mu_p*diff(p,x,2) + F*mu_p*diff(p,x) + r_m*p == G(x); %MARKED
P_ph = dsolve(eqn)
P_ph_int = int(P_ph,[350*10^-6 800*10^-6]) % need to integrate the P_ph in the interval of wavelength
P = P_dark+P_ph;
cond = [p(0) == N_v * exp(-phi_3/V_T), p(d_j) == N_v * exp(-phi_4/V_T)];
%This is the problem. I know how to handle it if the boundary is for P_ph or P_dark only. But right now, the boundary condition is used for P, the addition of the P_dark and P_ph. I can get the solution with undetermined constants but I can not use the boundary condition to find such constants.
P = dsolve(P,cond)
the bug info is
No differential equations found. Specify differential equations by using symbolic functions.
T = feval(symengine,'symobj::dsolve',sys,x,options);
sol = mupadDsolve(args, options);
0 Comments
Accepted Answer
Walter Roberson
on 16 Jul 2021
syms V_T mu_p F r_m % constant of ODE1
syms N_v phi_3 phi_4 d_j % constants for boundary condition
syms p(x) G(x) % p means hole density G(x) means generation rate
eqn = -V_T*mu_p*diff(p,x,2) + F*mu_p*diff(p,x) + r_m*p == 0; % equation in dark
P_dark(x) = dsolve(eqn)
eqn = -V_T*mu_p*diff(p,x,2) + F*mu_p*diff(p,x) + r_m*p -G(x) == 0; % equation in light (your suggestion) %MARKED
P_ph(x) = dsolve(eqn)
new_vars = setdiff(symvar(P_ph(x)), symvar(eqn))
P_total(x) = P_dark+P_ph % Solutions which boundary condition applied
eqn1 = P_total(0) == N_v * exp(-phi_3/V_T) %% boudary condition of P at x=0
eqn2 = P_total(d_j) == N_v * exp(-phi_4/V_T) %% boudary condition of P at x=d_j
eqns = [eqn1;eqn2]
symvar(eqns)
S = solve(eqns, new_vars)
temp = struct2cell(S);
new_vars(:) == vertcat(temp{:})
More Answers (1)
Torsten
on 16 Jul 2021
The equation for P reads
-V_T*mu_p*P''+ F*mu_p*P'+ r_m*P - G(x) = 0.
Now you can apply the boundary conditions.
3 Comments
Torsten
on 16 Jul 2021
In your code from above, P = P_dark + P_ph, not P = P_dark + P_ph_int.
What do you get as error message if you try to solve "my" differential equation for P with your boundary conditions ?
See Also
Categories
Find more on Equation Solving in Help Center and File Exchange
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!