In this example we use the package IntroToSymmetry.m to work out the point group of the equation 

Yxxx-f[x]*Yx^2=0

Clear all symbols in the current context.

ClearAll[Evaluate[Context[]<>"*"]]

First read in the package which is located in User Home Folder/Library/Mathematica/Applications/SymmetryAnalysis.

Needs["SymmetryAnalysis`IntroToSymmetry`"]

Enter the input equation as a string. Don't include the ==0 at the end.

inputequation=
"y'''[x]-f[x]*y'[x]^2";

The function y[x] 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.

rulesarray=
{"y'''[x]->f[x]*y'[x]^2"};

Enter the list of independent variables.

independentvariables={"x"};

Enter the list of dependent variables.

dependentvariables={"y"};

Enter the list of function names and constant names that need to be preserved when the equation is changed to (x1, y1) variables. In this case there are no such quantities and so we enter a null list.

frozennames={"f"};

Enter the maximum derivative order.

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-Bäcklund 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-Bäcklund groups and this can be useful when looking for contact symmetries.

xseon=1;

When searching for Lie-Bäcklund groups it is necessary to differentiate the input equation r times 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.

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

equationgenericvariables

invarconditiontable

The program expresses the determining equations in terms of z-variables. Here is the correspondence between z-variables and conventional variables.

ztableofrules

Here are the determining equations expressed in terms of z-variables. This output may be quite long. The equations in the table can be distinguished by the ==0 at the end of each item.

zdeterminingequations

Lets look at the determining equations in a slightly more readable form.

zdeterminingequationscolumn=TraditionalForm[Column[ToExpression[zdeterminingequations],{Left,Below}]]

Now solve the determining equations in terms of multivariable polynomials of a selected order. 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=3]]

Here are the infinitesimals for the independent variables expressed in terms of z-variables.

xsefunctions

and the infinitesimals for the dependent variables.

etafunctions

Express the xse functions in terms of x,y variables.

ToExpression[xsefunctions]/.ztableofrules

Express the eta function in terms of these variables.

ToExpression[etafunctions]/.ztableofrules

Here is a list of the various infinitesimal groups.

infinitesimalgroupsxy=infinitesimalgroups/.{z1->x,z2->y}

Column[infinitesimalgroupsxy]

The final result indicates that the equation is invariant under two groups. The parameter b10 corresponds to a translation in y. The parameter b14 needs further investigation.

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

First let's check to see if the b14 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 that generates the infinitesimal transformation of the pth derivative.

ToExpression[xsefunctions]/.{z1->x1,z2->y1[x1]}

ToExpression[etafunctions]/.{z1->x1,z2->y1[x1]}

n=Length[independentvariables]

m=Length[dependentvariables]

plocal is order of the desired derivative.

xse1noderiv=-b14 f[x1]/D[f[x1],x1]

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

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

invariancecondcoefficients={xse1noderiv,0,eta1x,0,eta1xxx}

invarconditiontable

invarianceconditionsum=
Expand[
invarconditiontable.invariancecondcoefficients]/.
{D[y1[x1],x1,x1,x1]->f[x1]*D[y1[x1],x1]^2}

Simplify[%]

Factor[%]

This final expression is equal to zero if and only if f[x] satisfies f[x]f''=(f')^2. Thus invariance of the equation holds if f=a*Exp[b*x]. Compare this with the exponential distribution of free stream velocity in the boundary layer over a flat plate.

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

Let's check.

f[x1_]:=a*Exp[b*x1]

ToExpression[xsefunctions]/.{z1->x1,z2->y1[x1]}

ToExpression[etafunctions]/.{z1->x1,z2->y1[x1]}

n=Length[independentvariables]

m=Length[dependentvariables]

xse1noderiv=-b14 f[x1]/D[f[x1],x1]

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

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

invariancecondcoefficients={xse1noderiv,0,eta1x,0,eta1xxx}

invarconditiontable

invarianceconditionsum=
Expand[
invarconditiontable.invariancecondcoefficients]/.
{D[y1[x1],x1,x1,x1]->f[x1]*D[y1[x1],x1]^2}

MaxMemoryUsed[ ]