BrianChung10 / SOP

2 stars 1 forks source link

Problems about solving the Bratu equation #5

Closed fz920 closed 2 years ago

fz920 commented 2 years ago

Hi, I am trying to plot the bifurcation diagram of the Bratu equation, but I get the SingularException Error for some values of lambda and now I get the MethodError: no method matching (::var"#F#17"{Float64})(::typeof(F)). The codes are in the Finite_Difference_1D / Goal_3_5.jl

Could you help with this? Thanks!

fz920 commented 2 years ago

MethodError: no method matching (::var"#F#17"{Float64})(::typeof(F)) Closest candidates are: (::var"#F#17")(::AbstractVector{T}) where T at ~/Desktop/group_project_deflation/SOP/Finite_Differences_1D/Goal_3-5.jl:89 Stacktrace: [1] newton(f::var"#F#17"{Float64}, x0::Function, max_iter::Int64, eps::Float64) @ Main ~/Desktop/group_project_deflation/SOP/Finite_Differences_1D/Goal_3-5.jl:35 [2] newton(f::Function, x0::Function) @ Main ~/Desktop/group_project_deflation/SOP/Finite_Differences_1D/Goal_3-5.jl:33 [3] bratu_solve(lambda::Float64, x0::Function) @ Main ~/Desktop/group_project_deflation/SOP/Finite_Differences_1D/Goal_3-5.jl:98 [4] top-level scope @ ~/Desktop/group_project_deflation/SOP/Finite_Differences_1D/Goal_3-5.jl:107

ioannisPApapadopoulos commented 2 years ago

Ignore my previous comment if you saw it. Tag me in issues like @fz920 otherwise I do not receive a notification. You have an error on line 104 n=length(lambda) but lambda is not defined. Just looking at the code you have sol = bratu_solve(lambda, F) in line 108 but have have the definition as function bratu_solve(lambda, x0) in line 87, i.e. you give it x0 in your definition but give it F in your call.

fz920 commented 2 years ago

@ioannisPApapadopoulos Thanks for pointing that out. I have fixed that error. But I found that for some particular values of lambda, I get different errors. For example, in the Finite_Difference_1D / Goal3_6.jl, when lambda = 1.1 I get the SingularException Error ,and when lambda = 2, I get the StackOverflow Error. Could you explain a bit about these errors? When ploting the bifurcation diagram, should we just ignore these values of lambda that cause the errors? Thanks in advance!

ioannisPApapadopoulos commented 2 years ago

SingularExceptionError means that the matrix you are inverting is numerically singular, i.e. not invertible. This is a nonlinear problem, so given some iterate x^k, you can imagine that the induced Jacobian could be singular in some circumstances, Newton is not perfect. You can try damping the Newton iteration (multiple the update with a constant <1) and see if that works.

However when I ran the code for lambda=1.1, I did not get an error at all, it worked?

Your StackOverflowError seems to come from the fact that deflated_newton command does not converge when lambda=2. Maybe try and damp that one?

If you cannot get values to work that's fine, just skip them. But ideally have enough values so that the bifurcation curve looks smooth.

dlfivefifty commented 2 years ago

If you do qr(A) \b it will do a least squares fit instead of trying to inverse… may be more robust

fz920 commented 2 years ago

Thank you for your suggestions. I have created a bifurcation diagram but the lines in the diagram are not smooth enough. I think some values in the sol are NaN caused this problem. Do you know how to improve this diagram? Thanks!

Screen Shot 2022-06-02 at 17 22 03
GianfrancoVieri commented 2 years ago

@ioannisPApapadopoulos @dlfivefifty , hi, may I know what is a global (or globalized) newton method? The example of the Painleve equation given in this https://epubs.siam.org/doi/pdf/10.1137/140984798 pdf mention about using a globalized newton algorithm to produce their graph which does converge and we have tried different values for the last row of our matrix but most of them do not converge. Thank you!

ioannisPApapadopoulos commented 2 years ago

Newton is a not globally convergent. Given a random initial guess, you're not guaranteed to converge to anything. If you want to make it global, then you also add a linesearch, i.e. an algorithm that adaptively damps the Newton update. Have you tried damping the Newton update with just a constant factor like 0.1?

GianfrancoVieri commented 2 years ago

Yeah, we tried it and we get infinity and -infinity for our values

GianfrancoVieri commented 2 years ago

yeah sorry, I asked in the wrong section, this should be under the painleve equation instead of the bratu equation

ioannisPApapadopoulos commented 2 years ago

This still applies to the disconnected Bratu bifurcation diagram you have. You also have a free choice of initial guess. Instead of just choosing a zero initial guess every time, you could use the solution(s) at the previous value of lambda as the guess for the next value of lambda. Much more likely to converge

dlfivefifty commented 2 years ago

This is called homotopy (where you vary a parameter and use the previous solution as an initial guess to the new problem).

@ioannisPApapadopoulos : is there a simple example of a discontinuous bifurcation example?

dlfivefifty commented 2 years ago

Newton is a not globally convergent

I’m pretty sure it’s globally convergent for convex functions. It would be interesting to know when it’s globally convergent when combined with deflation.

ioannisPApapadopoulos commented 2 years ago

Carrier's problem is an example of an ODE with a rich (and disconnected) bifurcation structure.

https://epubs.siam.org/doi/abs/10.1137/16M1096074?journalCode=smjmap

Sure! I mean "in general" it is not globally convergent. In particular all the problems we are looking at are nonconvex and so we might expect Newton to diverge from time to time.

You can construct more sophisticated Newton methods using linesearches and trust regions that are globally convergent for most problems.

Casper Beentjes' MSc thesis spent quite a bit of time looking at the convergence of Newton and how it might tie to deflation. https://cbeentjes.github.io/theses/2015-09-01-MSc-thesis

GianfrancoVieri commented 2 years ago

Oh I see, thank you so much! @ioannisPApapadopoulos @dlfivefifty , I'll try to reconfigure the code to adapt to the value of the last lambda.