sand15

797 Reputation

11 Badges

10 years, 10 days

MaplePrimes Activity


These are replies submitted by sand15

@Pepini

Look my reply to @Carl Love

@Carl Love

The figure at the end of the question cannot be obtained by rotating the red curve around some axis..

I don't know what is the name of this toroidal figure in English, discovered by the French mathematician Vincent Borrelli (Lyon University, 2012), but it is French dubbed the "Tore plat" (Flat Torus). This "Flat Torus" is an example of smooth fractals.
You can look here
https://www-sop.inria.fr/geometrica/events/wocg14-symmetry_and_periodicity/slides/thibert.pdf

https://www.google.com/search?q=square+flat+torus&tbm=isch&ved=2ahUKEwjQnYnH2q6AAxW9mCcCHcJFBAAQ2-cCegQIABAA&oq=square+flat+torus&gs_lcp=CgNpbWcQAzoFCAAQgAQ6CAgAELEDEIMBOgQIABADOggIABCABBCxAzoHCAAQigUQQzoLCAAQgAQQsQMQgwE6BAgAEB46BwgAEBMQgAQ6CAgAEAUQHhATOggIABAIEB4QE1CRFFiOPWDOQ2gAcAB4AIABcYgB8AiSAQQxNi4ymAEAoAEBqgELZ3dzLXdpei1pbWfAAQE&sclient=img&ei=vErCZNDPHL2xnsEPwosR&bih=808&biw=1729&client=firefox-b-d


This "Flat Torus" is related to the Nash–Kuiper theorem:
https://en.wikipedia.org/wiki/Nash_embedding_theorems

@Rana47 

  1. The ode has 3 parameters b, n and k.
     
  2. Its HPM counterparty depends on one more parameter wich is the order of the expansion (what I named the ansatz).

restart:

# In the sequel [REF] refers to the paper in your question

with(ColorTools):

#PDEtools[declare](f(x), prime = x);
#PDEtools[declare](v(x), prime = x);

ansatz := N -> g(x) = add((q^(i))*w[i](x), i = 0 .. N);

proc (N) options operator, arrow; g(x) = add(q^i*w[i](x), i = 0 .. N) end proc

(1)

HPMEq0 := n -> (1-q)*diff(g(x), x$2)+q*(diff(g(x), x$2)+n*b*(diff(g(x), x))^(n-1)*(diff(g(x), x$2))-k)

proc (n) options operator, arrow; (1-q)*(diff(g(x), `$`(x, 2)))+q*(diff(g(x), `$`(x, 2))+n*b*(diff(g(x), x))^(n-1)*(diff(g(x), `$`(x, 2)))-k) end proc

(2)

ansatz_order    := 3;
exponent_degree := 5;

# The initial guess corresponding to relation [REF] eq (29).

v[0]  := k/2*(x^2-2*x)+1;
Lv[0] := diff(v[0], x$2);

3

 

5

 

(1/2)*k*(x^2-2*x)+1

 

k

(3)

# Plug the ansatz into HPMEq0:

HPMEq1 := collect(eval(HPMEq0(exponent_degree), ansatz(ansatz_order)), q):

# Derive the successive edos according to [REF] eq (23), (25), (27) and
# generalize.

equ[0] := coeff(HPMEq1, q, 0) - diff(v[0], x$2) = 0;  # [REF] eq 23
for i from 1 to ansatz_order do
  equ[i] := coeff(HPMEq1, q, i) + diff(v[0], x$(i+1)) = 0
end do;

diff(diff(w[0](x), x), x)-k = 0

 

5*b*(diff(w[0](x), x))^4*(diff(diff(w[0](x), x), x))+diff(diff(w[1](x), x), x) = 0

 

20*b*(diff(w[0](x), x))^3*(diff(w[1](x), x))*(diff(diff(w[0](x), x), x))+5*b*(diff(w[0](x), x))^4*(diff(diff(w[1](x), x), x))+diff(diff(w[2](x), x), x) = 0

 

diff(diff(w[3](x), x), x)+5*b*(diff(w[0](x), x))^4*(diff(diff(w[2](x), x), x))+20*b*(diff(w[0](x), x))^3*(diff(w[1](x), x))*(diff(diff(w[1](x), x), x))+5*b*(2*(diff(w[0](x), x))^2*(2*(diff(w[0](x), x))*(diff(w[2](x), x))+(diff(w[1](x), x))^2)+4*(diff(w[0](x), x))^2*(diff(w[1](x), x))^2)*(diff(diff(w[0](x), x), x)) = 0

(4)

# Remark: equ[1] writes

coeff(HPMEq1, q, 1) + 'diff(v[0], x$(1+1))';

# and thus

value(%);

5*b*(diff(w[0](x), x))^4*(diff(diff(w[0](x), x), x))+diff(diff(w[1](x), x), x)-k+diff(v[0], `$`(x, 2))

 

5*b*(diff(w[0](x), x))^4*(diff(diff(w[0](x), x), x))+diff(diff(w[1](x), x), x)

(5)

# BCs given by eq (21) in the paper you present.

cond[0] := w[0](0) = 1, (D(w[0]))(1) = 0;

# BCs on the successive corrections

for j from 1 to ansatz_order do
  cond[j] := w[j](0) = 0, (D(w[j]))(1) = 0
end do;

w[0](0) = 1, (D(w[0]))(1) = 0

 

w[1](0) = 0, (D(w[1]))(1) = 0

 

w[2](0) = 0, (D(w[2]))(1) = 0

 

w[3](0) = 0, (D(w[3]))(1) = 0

(6)

answers := {dsolve({cond[0], equ[0]})};

for j from 1 to ansatz_order do
  ans := dsolve(eval({cond[j], equ[j]}, answers)):
  answers := answers union {ans}
end do:

{w[0](x) = (1/2)*k*x^2-k*x+1}

(7)

b_values := [seq(0.1..0.9, 0.2)]:
b_n      := numelems(b_values):
b_colors := EvenSpread(Color("RGB", "Red"), b_n):

S2 := eval(eval(rhs(ansatz(2)), answers), [q=1, k=1]);

plots:-display(
  seq(
    plot(
      eval(S2, b=b_values[b_i])
      , x=0..1
      , color=Color(b_colors[b_i])
      , legend=typeset('b'=b_values[b_i])
    )
    , b_i = 1..b_n
  )
  , title=typeset(["HPM solution", `<|>`(n=5, k=1, 'ansatz_order'=ansatz_order)])
  , view=[default, 0..1]
  , gridlines=true
);

(1/2)*x^2-x+1-(1/6)*b*(x-1)^6+(1/6)*b+(1/2)*b^2*(x-1)^10-(1/2)*b^2

 

 

# Numericl solution

SweepNumSol := proc(b_values, K, N)
  local b_n, b_colors, b_i, sys, sol, graphs:

  uses ColorTools, plots:

  b_n      := numelems(b_values):
  b_colors := EvenSpread(Color("RGB", "Red"), b_n):

  graphs := NULL:
  for b_i from 1 to b_n do
    sys    := eval({HPMEq0(N), g(0)=1, D(g)(1)=0}, [b=b_values[b_i], k=K, q=1]):
    sol    := dsolve(sys, numeric);
    graphs := graphs,
              odeplot(
                sol
                , [x, g(x)]
                , x=0..1
                , color=Color(b_colors[b_i])
                , legend=typeset('b'=b_values[b_i])
              )
  end do:
  display(
    graphs
    , title=typeset(["Numerical solution",`<|>`(n=5, k=1, 'ansatz_order'=ansatz_order)])
    , view=[default, 0..1]
    , gridlines=true
  )
end proc:

SweepNumSol(b_values, 1, 5)

 

 

Download parameterized_solution.mw

@Rana47 

You write "I have compared the ist and second orders problem but its look not same as reported in paper"
From the paper (equation 33)

  • green = zeroth order
  • red first order
  • blue second order

What I get

  • zeroth orders are the same
  • first orders are the same given n=5
  • second orders differ by a factor n in the paper.
    I noticed this in my mw file sayink that it could be either an error on my side or a typo in the paper.
    Here is the ode for the second order problem from paper ,eq 30.
    Here is this same ode same second in my file is  (remember n=5 and w <--> f, L(w2)=diff(w2, x$2))
    it is clear these two odes are identical.
    Given the doundary conditions I use are those of the paper (you can check that easily), and unless Maple is not capable to compute correctly the solution of this problem, I lean towards an error in the paper: the coefficient "n" in the second order term in equation 33 is a typo.

    In fact I am even sure there is an error in eq (33)  because after some substitutions the second order term is written this way (eq (34))

    As you can see the lultiplicative constant n in the second order term has disappeared.

This should close your first remark.

I dont understand what you say after "Sir, if we suppose ..."
Do you mean that "my" first order problem is not the one given by equations (25)-(26) ???
D
oes the absence of term k seem like a mistake on my side?
This is normal, this term is balanced by diff(v[0](x), x$2) which is equal to k.

Here is the same file as previously but with notations w instead of f and q instead of p to get closer to to make comparisons easier.
Comparison_made_easier.mw
Read it carefully before saying I did something wrong (which is of course perfectly possible)

@Rana47 

I've just seen that  @SHIVAS question has been removed and the the link
https://www.mapleprimes.com/questions/236755-How-To-Run-Or-Implement-HPM-Method-In-Maple

I sent you is obsolete.
Here is the file I attached to my answer to @SHIVAS

Download hpm_error_sand15.mw

@Rana47 

If you "do not know the other guys", why does your code contain the same names of variables, the same useless declare statement, the same way to write the ansatz, the same errors ?

"Please help to produce results mentioned in the pictures":

I do not inderstand what you are saying.
I answered
@SHIVAS a few hours ago and uploaded a code.
Go and see the answer here

https://www.mapleprimes.com/questions/236755-How-To-Run-Or-Implement-HPM-Method-In-Maple

Another advice: search for the name homotopy on this cite, you will find several answers.
To cite few
https://www.mapleprimes.com/questions/223044-Homotopy-Perturbation-Method
https://www.mapleprimes.com/questions/227515-Homotopy-Perturbation-Method
https://www.mapleprimes.com/questions/232586-How-To-Solve-System-Of-Nonlinear-ODE

@Rana47 

Did the three of you guys @Rana47 , @SHIVAS and @Madhukesh J K work together to ask the same question, or are you all the same person?

In any case I think that a correct attitude would be to reply the answers instead of opening a new thread on the same subject , with exactly the same bugged code.

Look here
https://www.mapleprimes.com/questions/236758-HPM-Error-Invalid-Input
and here
https://www.mapleprimes.com/questions/236747-Error-In-The-Program-HPM


As your code is strangely close to a recent code from another OP (are you and  @Madhukesh J Kthe same person?), I advice you to look here
https://www.mapleprimes.com/questions/236747-Error-In-The-Program-HPM
to see how HPM can be done.

@arashghgood 


You didn't answer all my remarks, for instance remark 2: where does a term like

(3/4)*k*omega*eta(k,t)*eta(k,t)

comes from?
It can come only from the Fourier transform of zeta(x, t) by its derivative wrt x:

FourierTransform(Int(zeta(x-z, t)*diff(zeta(z, t), z), z=-infinity..+infinity)) = I*k*eta(k, t)^2

FourierTransform(Int(zeta(x-z, t)*(diff(zeta(z, t), z)), z = -infinity .. infinity)) = I*k*eta(k, t)^2

(1)

 


Download Remark_2b.mw


Does the"main equation in real space" contain such a convolution product?

More suspect is the last term in your transformed ode: what is it term in real space where it come from?

diff(eta(k,t),t,t) + gamma*diff(eta(k,t),t) + omega^2*eta(k,t) + (3/4)*k*omega*eta(k,t)*eta(k,t) + (3/2)*k*eta(k,t)*(omega*eta(k,t) + (3/4)*k*eta(k,t)*eta(k,t)) + (omega*gamma*eta(k,t) + (3/2)*k*eta(k,t)*(gamma*eta(k,t)))*(diff(eta(k,t),t) + gamma*eta(k,t))/(omega*eta(k,t) + (3/4)*k*eta(k,t)*eta(k,t)) = 0:
op(-1, lhs(%))

(omega*gamma*eta(k, t)+(3/2)*k*eta(k, t)^2*gamma)*(diff(eta(k, t), t)+gamma*eta(k, t))/(omega*eta(k, t)+(3/4)*k*eta(k, t)^2)

(1)

 


Download Remark_3.mw

I think you should provide us the "main equation in real space" with its initial and/or boundary conditions.

in the system you dsolve.
Did you mean equa?

@acer 

@Hullzie16 code worked perfect well with Maple 2015.2.
Maybe this information will help the development team to find where the bug comes from.

@dharr 

I voted up, but I believe there are a few points you could consider.

Firstly I don't thik writting K[3] = K[1]*k[3] is appropriate because neither K[1]nor K[3] intervene by themselves in the equation:  eqn only depends on the difference K[3]-K[1].
I think it would have been better to set K[3] = K[1] + k instead, which reduces the number ot parameters by 1.

Next, dimension analysis shows that gamma[1] and gamma[2] have the same dimension: this could lead to a simpler final equation with still one less parameter.

(last point: it is easy to show that theta is dimensionless).


Completing @acer's reply:

  • You meant probably k*alpha instead of ?
    Copy-pasting into a 1D mode worksheet mode
    `k&alpha;`
  • There is no need do declare GAMMA as local as soon as you define GAMMA... like the original (when its argument is an integer).
    By the way
    GAMMA(k*alpha+1) / GAMMA(k*alpha+alpha+1) = GAMMA(k+1) / GAMMA(k+2) = 1/(k+1)
    
  • Within your solve command, the first argument is
    1/(k+1)*(diff(U[k](x), x, x)+2*sum(U[r](x)*(diff(U[k-r](x), x)), r = 0 .. k)-(diff(sum(U[r](x)*U[k-r](x), r = 0 .. k), x)))
    

    (
    I used sum instead of add to "show" what this term [TERM] looks like

    )
    What result do you expect by doing

    solve(TERM, U[k+1](x))

    given TERM doesn't contain U[k+1](x) ?

@raj2018 

Integrating formally f1^4 is a dead end for the same reason that integrating f2 is.
f1^4 contains terms of the form:

(1-x/a)^(-kc+1/2)*(1-x/b)^(-kh+1/2)

So the only possibility is to use a numerical integration.



Several errors:

  1. No termination of the loop (end do is missing).
  2. What do you expectto do writting f[i] = f[i+1] ?
  3. n being undefined there is no chance that an hypothetical result matches the expected one (which doesn't contain n).
    Are those n typos?
     

Point 1 put apart, your problem has nothing to do with Maple: this is more a question about how you are capable or not to write correctly your recurrence relations.
Fix all these points and come again with something more consistent if you still get errors.

First 6 7 8 9 10 11 12 Last Page 8 of 25