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]> 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 :
>
>
>
>
