Ok, I understand. It works and it is really awesome. Thank you !
Le samedi 21 juin 2014 12:17:22 UTC+2, Sheehan Olver a écrit : > > Hi, > > > You didn’t quite get the Newton iteration right: if you want to solve > > B u + [1,0] = 0 > L u + g(u) = 0 > > then Newton iteration becomes > > u = u - [B, L + g’(u)]\[B u + [1,0], L u + g(u)] > > > i.e., your bc right hand side is not right. Below is the corrected code. > > Cheers, > > Sheehan > > > > x=Fun(identity, Interval(0.,1.)) > d=x.domain > B=neumann(d) > D=diff(d) > # Solves Lu + g(u)-1==0 > L = D^2 + 2/(x.+1)*D > g = u -> -(exp(u)-exp(-u)); gp = u -> -(exp(u)+exp(-u)) > > u=-0.3*x+0.5 #initial guess > > for k=1:5 # this crashes if Ii increase the number of iterations > u=u-[B,L+gp(u)]\[diff(u)[0.]+1.,diff(u)[1.],L*u+g(u)]; > end > > plot(u) > > On 21 Jun 2014, at 6:54 pm, 'Stéphane Laurent' via julia-users < > [email protected] <javascript:>> wrote: > > Hello, > > I'd like to solve this equation with Neumann boundary conditions. My code > below does not work. Am I doing something bad or is it a failure of the > Newton algorithm ? > > *# solves u" = (exp(u)-exp(-u))-2u'/(x+1) , u'(0)=-1, u(1)=0* > *x=Fun(identity, Interval(0.,1.))* > *d=x.domain* > *B=neumann(d)* > *D=diff(d)* > *# Solves Lu + g(u)-1==0* > *L = D^2 + 2/(x.+1)*D* > *g = u -> -(exp(u)-exp(-u)); gp = u -> -(exp(u)+exp(-u))* > > *u=-0.3*x+0.5** #initial guess * > > *for k=1:5 # this crashes if Ii increase the number of iterations* > * u=u-[B,L+gp(u)]\[-1.,0.,L*u+g(u)];* > *end* > > > The solution should look like that : > > > <https://lh5.googleusercontent.com/-BCYNzh8ckCA/U6VIJcT2-tI/AAAAAAAAAGI/r08hknI2EIg/s1600/Screenshot+from+2014-06-21+10%3A37%3A46.png> > > >
