Solving a system of differential equations in Matlab

Question

Consider following system of differential equations:

$$\begin{aligned} X' &= Z \cdot (2.5 - |Y| + \cos(\Phi))\\ Y' &= Z \cdot(\sin(\Phi))\\ \Phi' &= Z \cdot(\sin^2(\Phi) + \mbox{sgn}(Y) \cdot \cos^2(\Phi))\\ Z' &= 0\end{aligned}$$

I have the boundary conditions

$$X(0)=x_0, \qquad Y(0)=y_0, \qquad X(1)=x_1, \qquad Y(1)=y_1$$

1. How do I solve this in Matlab numerically?

2. How do I solve this analytically?

My attempt:

xmesh = linspace(0,1,1000);
solinit = bvpinit(xmesh,@guess);

sol = bvp4c(@bvpfun, @bcfun, solinit);
k=figure(2)
plot(sol.x,sol.y,'-o')

function dydx = bvpfun(t,y) 
dydx = [y(4,:)*(2.5-abs(y(1,:))+cos(y(3,:)))
       y(4,:)*(0+sin(y(3,:)))
       y(4,:)*(0*(sin(y(3,:)))^2-(-sign(y(2,:)))*(cos(y(3,:)))^2+(0)*sin(y(3,:))*cos(y(3,:)))
       0];
end

function res = bcfun(ya,yb)
res = [ya(1)-1
       ya(2)+2.5
       yb(1)-1
       yb(2)-2.5];
end

function  u = guess(t)
u = [1
     -2.5+t*5
     30
     15];
end

The Problems i have :

1. the solution depends strong on the guess and i want to avoid that.
2. the last variable $Z$ (constant) should be interpreted as a time, but in some of my computations it is negative.

is there a method for solving these equations with respect to the boundary conditions without a guess???

Another try with shooting method:

function shooting_method
clc
clear all
x=[2 1]
options=optimset('Display','iter')
x1=fsolve(@solver,x)
end

function L=solver(x)
options=odeset('RelTol', 1e-8,'AbsTol',[1e-8 1e-8 1e-8 1e-8]);
[t,u]=ode45(@equation,[0 1],[1 x(1) -2.5 x(2)], options);
s=length(t)
L=[u(s,1)-1, u(s,2)-2.5];
figure(3)
plot(t,u(1,:),u(2,:),u(3,:),u(4,:))
end

function dydt = equation(t,y)
dy=zeros(4,1)
dy(1) = y(4)*(2.5-abs(y(1))+cos(y(3)))
dy(2) = y(4)*(sin(y(3)))
dy(3) = y(4)*(-(-sign(y(2)))*(cos(y(3)))^2)
dy(4) = 0     
end

Here i get the error:

Error in shooting_method>solver (line 11)
[t,u]=ode45(@equation,[0 1],[1 x(1) -2.5 x(2)], options);

Error in fsolve (line 242)
            fuser = feval(funfcn{3},x,varargin{:});

Error in shooting_method (line 6)
x1=fsolve(@solver,x)

Caused by:
    Failure in initial objective function evaluation. FSOLVE cannot continue.

Answer 1

You can also use the "unknown parameters" mechanism dydx = bvpfun(t,y,p), res = bcfun(ya,yb), solinit = bvpinit(xmesh,@guess,[Z0]);, to avoid adding the additional component for $Z$.

A non-linear boundary value problem can be expected to have an arbitrary number of solutions, from zero to infinitely many. You will need to have a good idea of the shape of the desired solution to pass to bvpinit to fix on the exact solution close to it.

The standard in BVP solvers is to use collocation methods or at least implicit methods of a high order. These react very badly to discontinuities as introduced with the signum of $Y$. Even the absolute value might lead to an increased mesh density around points where it is zero, which might lead to exceeding the allowed number of mesh nodes. It is worth a try to replace them with smooth mollifications.