lam3 = 1e-6;
deff = 20e-12;
L=2e-2;
n=2.2;
c = 299792458;
e0 = 8.854187e-12;

phispan = [0 2*pi];

Ip = 100e7;
Is = 10e-2;

ais = sqrt((2*Is)/(n*e0*c));
aip = sqrt((2*Ip)/(n*e0*c));


aInit = [ais aip];

[P,A] = ode45(@cmeq,phispan,aInit);



function dA = cmeq(phi,A)

lam3 = 1e-6;
deff = 20e-12;
L=2e-2;
n=2.2;
c = 299792458;

w2 = c/lam3;
w1 = w2/2;

K1 = -(1i*w1*c*deff)/n;
K2 = -(1i*w2*c*deff)/n;


dA = zeros(2,1);

dA(1) = K1*A(2)*(A(1))';
dA(2) = K2*(A(1))^2;

A(1) = abs(A(1))*exp(1i*phi);

end