1
$\begingroup$

I need to set a "boundary" condition not at the boundaries of the computing domain but inside the domain during solving an ODE with FDM. The problem is a boundary value problem, which has been implemented based on pdetoae devised by @xzczd.

  1. Prepare the mesh:

    L = 6; domain = {-L, L};
    points = 21;
    difforder = 4; grid = Array[# &, points, domain];
    (*please find the function from the link above*)
    ptoafunc = pdetoae[u[x], grid, difforder];
    
  2. Discretize the system. Please note that the last condition u'[0] == 0 of bc is what I want to impose inside the domain [-L,L]:

    eq = u[x] u'[x] + u''[x] + u''''[x] - 2*u'[x] == 0;
    bc = {u[-L] == 0, u[L] == 0, u'[-L] == 0, u'[L] == 0, u'[0] == 0};
    
    del = #[[3 ;; -3]] &;
    ae = ptoafunc@eq // del;
    aebc = ptoafunc@bc;
    
  3. Solve the resulting system with FindRoot:

    ic[x_] = 10*Cos[\[Pi]/(2 L) x];
    solrule = FindRoot[{ae, aebc}, Table[{u@x, ic@x}, {x, grid}]]
    

Mathematica gives the following error:

FindRoot::nveq: "The number of equations does not match the number of variables in FindRoot[...]."

I also tried to remove one or two boundary conditions among the first four of bc, for example, u[-L] == 0,(*u[L] == 0,*)u'[-L] == 0,(*u'[L] == 0,*). However, the above error remains. Here, without the last condition u'[0] == 0, the code works well, but I need to additionally set this "boundary" condition to obtain a special function $u(x)$.

Thank you for any suggestion!

Update (according to xzczd's suggestion): Such a b.c could be imposed with

midcond = ptoafunc[u'[0] == 0][[(points + 1)/2]]

Now we would have 5 b.c.s, and the problem appears to be ill-posed, so the code gives the same warning again. Thus, we must remove one of them. Next, using aebc = ptoafunc@bc~Join~{midcond} as the new b.c.s. FindRoot gives a reasonable solution.

$\endgroup$
  • $\begingroup$ Currently pdetoae cannot handle b.c. that's inside the domain, it'll just discretize u'[0]==0 as if it's u'[x]==0, just execute ptoafunc[u'[0]==0] and observe. It's not hard to circumvent this, you can use Part to extract the needed b.c. from ptoafunc[u'[0]==0]. However, are you sure these 5 b.c.s form a well-posed problem? Notice imposing b.c. inside the domain is essentially equivalent to split the original domain. $\endgroup$ – xzczd Jun 18 at 6:54
  • $\begingroup$ @xzczd thank you for the help. Have a quick review for my update? $\endgroup$ – jsxs Jun 18 at 7:40
  • $\begingroup$ Removing one b.c. at $x=\pm L$ is a possible way to avoid the warning of course. Now there's nothing wrong in the programming side, but I think one needs to check carefully to see if the result is reasonable, because the b.c. imposed in this case is a bit unusual, I should say. $\endgroup$ – xzczd Jun 18 at 7:49
  • $\begingroup$ There is an obvious solution u[x]=0. $\endgroup$ – Alex Trounev Jun 18 at 11:38
  • $\begingroup$ @AlexTrounev yes, it admits the trivial solution u=0. However, I want to find the non-trivial ones :) $\endgroup$ – jsxs Jun 18 at 12:28

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Browse other questions tagged or ask your own question.