Solving for Coefficient in Nonlinear 2nd-Order Differential Equation

1 view (last 30 days)
Hi,
I'm trying to numerically solve for the coefficient "a" in a nonlinear second-order differential equation of the form:
y''(r) = y(r)(a - b*abs(y(r))/r + (c*r^3 - d)*y(r)^2)
given boundary conditions:
y(0) = y(inf) = 0
integral from 0 to inf of y(r)^2 dr = 1
However I don't seem to be able to find the tools to do this neatly. I have another platform that can do this, but I'd very much like to cross-check the results in Matlab!
Any help is appreciated!

Answers (2)

Torsten
Torsten on 5 Jan 2018
Edited: Torsten on 5 Jan 2018
Use "fzero" or "fsolve" to solve the equation
(integral_{0}^{inf} y(a,r)^2 dr) - 1 = 0
for a.
To this end, in the function routine from "fzero" or "fsolve", you will be given a suggestion for a. With this value for a, you can use "bvp4c" to solve the boundary value problem
y''(r) = y(r)(a - b*abs(y(r))/r + (c*r^3 - d)*y(r)^2)
given boundary conditions:
y(0) = y(inf) = 0
Once you have its solution in discrete points, use "trapz" to approximate integral_{0}^{inf} y(r)^2 dr and return integral_{0}^{inf} y(r)^2 dr - 1 to "fzero" or "fsolve".
I think one problem with a numerical solution will be that your boundary value problem always has y = 0 for all r as solution.
Best wishes
Torsten.
  1 Comment
Alexander Herrod
Alexander Herrod on 21 Jul 2018
Thanks for the reply, and sorry it's taken so long!
I don't assume the dependency of y(r) on any of the constants. The integration requirement is only to be applied at the end.
If I try your suggestion, Matlab throws me an error that it can't find and solutions to that boundary condition (because of course there are infinite possible functions!)

Sign in to comment.


John BG
John BG on 5 Jan 2018
Edited: John BG on 15 Jan 2018
Hi Alexander
1.
One of the tools available is MuPAD, start MuPAD with
mupad
once open define the equation, but avoid the abs(), it only complicates the result and since r seems to be a cylindrical coordinate, if r is the radius, then there is no point that r takes negative values
.
the result after ode:solve(o) is what comes up if including the abs()
the resulting functions do not seem that useful.
2.
It turns out that the Jacobi elliptic functions are solution to a similar differential equation:
clear all;clc;close all
rspan=[-20 20];
dr=0.0001;
r=[rspan(1):dr:rspan(2)];
M = 0.7;
[S,C,D] = ellipj(r,M);
figure(1);plot(r,S,r,C,r,D);
legend('SN','CN','DN','Location','best');grid on;
title('Jacobi Elliptic Functions sn,cn,dn');
a=1;b=1;c=.1;d=.1;
vS=d2ydr2(r,S,a,b,c,d);
vC=d2ydr2(r,C,a,b,c,d);
vD=d2ydr2(r,D,a,b,c,d);
d2Sdt2=diff(diff(S));
d2Cdt2=diff(diff(C));
d2Ddt2=diff(diff(D));
errS=d2Sdt2-vS([3:end]);
errC=d2Cdt2-vC([3:end]);
errD=d2Ddt2-vD([3:end]);
figure(2);plot(r([3:end]),abs(errS),r([3:end]),abs(errC),r([3:end]),abs(errD));grid on;
3.
now M=.99 CN and DN (not the blue one) show far smaller error
rspan=[1 5];
dr=0.0001;
r=[rspan(1):dr:rspan(2)];
M = 0.99;
[S,C,D] = ellipj(r,M);
vS=d2ydr2(r,S,a,b,c,d);
vC=d2ydr2(r,C,a,b,c,d);
vD=d2ydr2(r,D,a,b,c,d);
d2Sdt2=diff(diff(S));
d2Cdt2=diff(diff(C));
d2Ddt2=diff(diff(D));
errS=d2Sdt2-vS([3:end]);
errC=d2Cdt2-vC([3:end]);
errD=d2Ddt2-vD([3:end]);
figure(3);plot(r([3:end]),abs(errS),r([3:end]),abs(errC),r([3:end]),abs(errD));grid on;
if you find this answer useful would you please be so kind to consider marking my answer as Accepted Answer?
To any other reader, if you find this answer useful please consider clicking on the thumbs-up vote link
thanks in advance for time and attention
John BG
  1 Comment
Alexander Herrod
Alexander Herrod on 22 Jul 2018
Thanks John,
Sorry it's taken me so long to get to this. Is there any way to input boundary conditions of the kind I've given?
Thanks again!

Sign in to comment.

Community Treasure Hunt

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

Start Hunting!