How do I make this substitution in Mathematica?

Question

I'm just getting started with Mathematica, and I've got what must be a pretty basic question about making substitutions but I can't get it to work.

I'd like to find the Euler-Lagrange equations for a functional of the function phi[x,y] and then make a substitution for the function phi[x,y]

If I enter the following:

VariationalD[tau*phi[x, y]^2 - 2*phi[x, y]^4 + phi[x, y]^6 + Dot[D[phi[x, y], {{x, y}}], D[phi[x, y, {{x, y}}]]], phi[x, y], {x, y}]

I get

Plus[Times[2,tau,phi[x,y]],Times[-8,Power[phi[x,y],3]],Times[6,Power[phi[x,y],5]],Times[-2,Plus[Derivative[0,2][phi][x,y],Derivative[2,0][phi][x,y]]]]

Now if I try % /. phi[x,y] -> phi0[x,y] + psi[x,y] it makes the substitution for all the polynomial terms, but not for the Derivative terms.

How do I force the substitution into those functions?

+1, for making me struggle to come up with a simple, generic way to accomplish this. — rcollyer, Dec 14 '10 at 21:38
This may be obvious to physicists, but what substitution are you trying to do for `Derivative[0, 2][phi][x, y]` term? — Yaroslav Bulatov, Dec 14 '10 at 23:09
@Yaro phi[x,y]->phi0[x,y]+psi[x,y] ... probably he is doing some perturbation analysis. After that you keep only the most significant terms is psi. That's the way ODEs always give you an harmonic oscillator ... or chaos. That's physics :) — Dr. belisarius, Dec 15 '10 at 00:10

score 5 · Accepted Answer · edited May 23 '17 at 12:11

I agree with all of what rcollyer says, but I think his final solution might be a little opaque.

The simplest rule that I could come up with (which is the basically the same as rcollyer's) is

{phi[x__] :> phi0[x] + psi[x], f_[phi][x__] :> f[phi0][x] + f[psi][x]}

or something with less possible side effects is

{phi[x__] :> phi0[x] + psi[x], Derivative[n__][phi][x__] :> Derivative[n][phi0][x] + Derivative[n][psi][x]}

It would be a lot easier if Derivative had a Default property (compare Default[Times] with Default[Derivative]). It should be something like Default[Derivative] := Sequence[] but unfortunately that doesn't play nice with the pattern matching.

Getting back to your question, you probably want to define something like

VariationalD[expr_, sym_, var_] := Module[{
  vRule = {sym[x__] :> sym[x] + var[x], 
    Derivative[n__][sym][x__] :> Derivative[n][sym][x] + Derivative[n][var][x]}}, 
  (expr /. vRule) - expr]

Where the variation var of the symbol sym is assumed to be small. Of course what you then need to do is series expand around var=0 and only keep the linear part. Then use integration by parts on any term which has derivatives of var. All of which should be included in the above module.

+1, interesting answer. My answer is opaque, however I object to the need to create a special rule every time I need to substitute in the sum of a function. I'd rather have a rule I can use without modification, but to each their own. Your rule would be easier to craft off the top of my head. — rcollyer, Dec 15 '10 at 03:08
As to your point on the possible side effects of my answer, the pattern `a_[b__][c__]` will only match those heads where there are two or more sets of square brackets. Additionally, `Distribute[a[b]]` will only actually rearrange things if `b` contains `Plus`, and it can be no more than one level down, i.e. `f[a,x+y,b]` becomes `f[a,x,b]+f[a,y,b]` but `f[a,g[x+y],b]` is unaffected. `Through` is the only potential problem child in that `Through[a b[x]]` is unevaluated, so the pattern could be improved, slightly. — rcollyer, Dec 16 '10 at 22:23

score 2 · Answer 2 · edited May 23 '17 at 12:11

First, you misplaced a ] in your second derivative term, it should read D[phi[x, y], {{x, y}}]] not D[phi[x, y, {{x, y}}]]].

That said, replacement in Mathematica can be tricky, as has been pointed out in other questions. That isn't to say it is impossible, just requires some work. In this case, the problem comes in in that phi[x,y] is different from Derivative[2, 0][phi][x, y]. So, your pattern won't match the derivative term. The simplest thing to do is to add the rule

Derivative[a__][phi][x__]:> Derivative[a][phi0][x] + Derivative[a][psi][x]

to your list of replacement rules. Three things to note: 1) I use ReplaceDelayed so that both types of derivatives will match without writing multiple rules, 2) since I can use patterns, I named them so that I can refer to them on the RHS of the rule, and 3) I used a double underscore when defining a and x which will match one or more items in a sequence.

Of course, that isn't the most satisfying way to approach the problem, as it will require you to write two rules every time you wish to this sort of replacement. It turns out a more general approach is surprisingly difficult to accomplish, and I'll have to get back to you on it.

Edit: This requires a double replacement, as follows

<result> /. phi -> phi0 + psi /. a_[b__][c__] :> Through[Distribute[a[b]][c]]

Distribute ensures that the derivative works correctly with Plus, and Through does the same with the function args c. The key is that the Head of Derivative[2, 0][phi][x, y] is Derivative[2, 0][phi], hence the several levels of square brackets in the rule.

Thanks for the answer. It works, but I'm having some trouble still with parsing it. What exactly is the pattern a_[b__][c__] matching? I guess for example a_ could match Derivative, b__ could match [2,0], and c__ could match [phi0+psi]. So then the rhs of the RuleDelayed would be Through[Distribute[Derivative[2,0]][phi0+psi]]. I must be interpreting the brackets wrong or the order of operations wrong, but it looks to me like the Distribute function doesn't have anything to distribute Derivative over. Unless maybe the Through function is being evaluated first? — 2daaa, Dec 16 '10 at 01:42
@Ranjit, that was the hard one. This pattern, `a_[b__][c__]` is easier to read from right to left: `[c__]` matches anything with square brackets at the end, `[b__]` then matches the 2nd to last square brackets, so `a_` matches everything else. So, `a_ :> Derivative[2,0]`, `b__ :> phi0 + psi`, and `c__ :> Sequence[x,y]`. Then, the purpose of `Distribute` is to interchange `Derivative[2,0]` and `Plus`, and `Through` attaches `[x,y]` to each derivative. — rcollyer, Dec 16 '10 at 02:47

How do I make this substitution in Mathematica?

2 Answers2