In this example we use the package IntroToSymmetry.m to work out the point group of the boundary layer stream function equation, 

ΨyΨxy-ΨxΨyy-UeUex-nu*Ψyyy=0.

This next command turns off a bunch of spurious spelling error messages.

Off[General::spell]

First read in the package which is located in a folder called SymmetryAnalysis inside the Mathematica folder.

Needs["SymmetryAnalysis`IntroToSymmetry`"]

Enter the input equation as a string. Don't include the == 0 at the end. Notice we are treating Ue as a dependent variable, in effect making the system underdetermined.

inputequation=
"D[Ψ[x,y],y,y,y]+(1/nu)*Ue[x,y]*D[Ue[x,y],x]+(1/nu)*D[Ψ[x,y],x]*D[Ψ[x,y],y,y]-
(1/nu)*D[Ψ[x,y],y]*D[Ψ[x,y],x,y]";

The function Ψ[x,y] is a solution of the equation and this constraint must be applied in the form of a rule to the invariance condition. Be careful to check signs. Note the additional rule that makes Ue independent of y.

rulesarray=
{"D[Ψ[x,y],y,y,y]->-(1/nu)*Ue[x,y]*D[Ue[x,y],x]-(1/nu)*D[Ψ[x,y],x]*D[Ψ[x,y],y,y]+
(1/nu)*D[Ψ[x,y],y]*D[Ψ[x,y],x,y]",
"D[Ue[x,y],y]->0"};

Enter the list of independent variables.

independentvariables = {"x", "y"};

Enter the list of dependent variables.

dependentvariables={"\!\(\*
StyleBox[\"Ψ\",\nFontFamily->\"Symbol\",\nFontSize->14,\nFontWeight->\"Bold\"]\)","Ue"};

Enter the list of function and constant names that must be preserved when the equations are expressed in terms of generic variables.

frozennames = {"nu"};

Enter the maximum derivative order of the input equation(s).

p = 3;

The maximum derivative order that the infinitesimals are assumed to depend on is specified by the input parameter r. This parameter is only nonzero when the user is looking for Lie contact groups or Lie-Backlund groups. For the usual case where one is searching for point groups set r=0.

r = 0;

When searching for Lie-Backlund groups (r=1 or greater) one can, without loss of generality, leave the independent variables untransformed. The corresponding infinitesimals (the xse's) are set to zero by setting xseon=0. If one is searching for point groups then set xseon=1. The choice xseon=1 is also an option when looking for Lie-Backlund groups and this can be useful when looking for contact symmetries.

xseon = 1;

When searching for Lie-Backlund groups it is necessary to differentiate the input equation with respect to each of the independent variables producing derivatives of order p+r. These higher order differential consequences are appended to the set of rules applied to the invariance condition. This process is carried out automatically when internalrules=1. For point groups the equation or equation system is the only rule or set of rules needed and one sets internalrules=0.

internalrules = 0;

Now work out the determining equations of the Lie point group that leaves the equation invariant. The output is available as a table of strings called zdeterminingequations. The Mathematica function Timing outputs the time required for the FindDeterminingEquations function to execute.

Timing[
      FindDeterminingEquations[
independentvariables, 
        dependentvariables, frozennames, p, r, 
        xseon, 
inputequation, rulesarray, 
        internalrules]]

equationgenericvariables

invarconditiontable

zdeterminingequations;

Length[zdeterminingequations]

ztableofrules

Now solve the determining equations in terms of multivariable polynomials. The Mathematica function Solve uses Gaussian elimination to solve a large number of linear equations for the polynomial coefficients. The time roughly follows,

time/timeref=((number of equations)/(number of equationsref))^n

where the exponent is between 2.4 and 2.7. The Mathematica function Timing outputs the time required for the SolveDeterminingEquations function to execute.

Timing[SolveDeterminingEquations[independentvariables,dependentvariables,r,xseon,
zdeterminingequations,order=4]]

Here are the infinitesimal transformation functions for the independent variables.

xsefunctions

and for the dependent variables.

etafunctions

infinitesimalgroupsxy = 
      infinitesimalgroups /. 
{z1 -> x, z2 -> y, 
          z3 -> psi, z4 -> Ue}

Check the groups; xse's are on the left, eta's are on the right.

ColumnForm[infinitesimalgroupsxy]

Note the higher powers of x that appear in the transformation of the wall normal coordinate. This suggests the existence of an arbitrary function of x in the transformation of y.

infinitesimalgroupsxy1 = Drop[infinitesimalgroupsxy, {2, 6}]

MakeCommutatorTable[independentvariables, dependentvariables, 
      infinitesimalgroupsxy1]

MatrixForm[commutatortable]

*********************************************************

First let's check to see if these infinitesimals satisfy the invariance condition by forming all the terms that appear in the invariance condition. Here we will make use of the function PthInfinitesimal the generates the infinitesimal transformation of the pth derivative.

ToExpression[xsefunctions] /. {z1 -> x1, z2 -> x2, 
        z3 -> y1[x1, x2], z4 -> y2[x1, x2]}

ToExpression[etafunctions]/.
{z1->x1,z2->x2,z3->y1[x1,x2],z4->y2[x1,x2]}

n = Length[independentvariables]

m = Length[dependentvariables]

The parameter plocal is order of the desired derivative.

eta1x=PthInfinitesimal[
n,m,yindex=1,plocal=1,r,xseon,j1=1]

eta1y=PthInfinitesimal[
n,m,yindex=1,plocal=1,r,xseon,j1=2]

eta1xy=PthInfinitesimal[
n,m,yindex=1,plocal=2,r,xseon,j1=1,j2=2]

eta1yy=PthInfinitesimal[
n,m,yindex=1,plocal=2,r,xseon,j1=2,j2=2]

eta1yyy=PthInfinitesimal[
n,m,yindex=1,plocal=3,r,xseon,j1=2,j2=2,j3=2]

eta2noder = (-a25  +  b115)*y2[x1, x2]

eta2x = 
      PthInfinitesimal[
n, m, yindex = 2, plocal = 1, r, 
        xseon, j1 = 1]

invariancecondcoefficients={
0,0,0,eta2noder,eta1x,eta1y,0,eta1xy,eta1yy,0,0,0,eta1yyy,
eta2x,0,0,0,0,0,0,0,0}

invarianceconditionsum=
Expand[
invarconditiontable.invariancecondcoefficients]/.
{D[y1[x1,x2],x2,x2,x2]->
-(1/nu)*y2[x1,x2]*D[y2[x1,x2],x1]-
(1/nu)*D[y1[x1,x2],x1]*D[y1[x1,x2],x2,x2]+
(1/nu)*D[y1[x1,x2],x2]*D[y1[x1,x2],x1,x2],
D[y2[x1,x2],x2]->0}

Simplify[%]

***********************************************

The presence of a truncated power series with order equal to the order of the trial series in the expression for xse2 suggests that this infinitesimal may depend on an arbitrary function of x1. Lets try this out and see if the invariance condition is still satisfied. The arbitrary function is g[x]

xse1[x1_,x2_,y1[x1_,x2_],y2[x1_,x2_]]:=a10 + (a25 + b115)*x1

xse2[x1_,x2_,y1[x1_,x2_],y2[x1_,x2_]]:=g[x1]+ a25*x2

eta1[x1_,x2_,y1[x1_,x2_],y2[x1_,x2_]]:=b10 + b115*y1[x1,x2]

eta2[x1_,x2_,y1[x1_,x2_],y2[x1_,x2_]]:=(-a25 + b115)*y2[x1,x2]

eta2eval = (-a25 + b115)*y2[x1, x2]

etaa1x=PthInfinitesimal[
n,m,yindex=1,plocal=1,r,xseon,j1=1]

etaa1y=PthInfinitesimal[
n,m,yindex=1,plocal=1,r,xseon,j1=2]

etaa1xy=PthInfinitesimal[
n,m,yindex=1,plocal=2,r,xseon,j1=1,j2=2]

etaa1yy=PthInfinitesimal[
n,m,yindex=1,plocal=2,r,xseon,j1=2,j2=2]

etaa1yyy=PthInfinitesimal[
n,m,yindex=1,plocal=3,r,xseon,j1=2,j2=2,j3=2]

etaa2x=PthInfinitesimal[
n,m,yindex=2,plocal=1,r,xseon,j1=1]

invariancecondcoefficientsnew={
0,0,0,eta2eval,etaa1x,etaa1y,0,etaa1xy,etaa1yy,0,0,0,etaa1yyy,
etaa2x,0,0,0,0,0,0,0,0}

invarianceconditionsumnew=
Expand[
invarconditiontable.invariancecondcoefficientsnew]/.
{D[y1[x1,x2],x2,x2,x2]->
-(1/nu)*y2[x1,x2]*D[y2[x1,x2],x1]-
(1/nu)*D[y1[x1,x2],x1]*D[y1[x1,x2],x2,x2]+
(1/nu)*D[y1[x1,x2],x2]*D[y1[x1,x2],x1,x2],
D[y2[x1,x2],x2]->0}

Simplify[%]

So the general boundary layer streamfunction equation admits four simple point groups,

a10:  (xse1,xse2,eta1,eta2)=(1,0,0,0)
b10:  (xse1,xse2,eta1,eta2)=(0,0,1,0)
a25:  (xse1,xse2,eta1,eta2)=(x,y,0,-Ue)
b115: (xse1,xse2,eta1,eta2)=(x,0,Ψ,Ue)

In addition the equations admit the infinite dimensional translation group (xse1,xse2,eta1,eta2)=(0,g[x],0,0).

MaxMemoryUsed[ ]