In this example we use the package IntroToSymmetry.m to work out nonclassical symmetries of the following form of the Fokker-Planck equation, 

Uxx+xUx-Ut=0.

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.

\(\(fpequation = 
"D[u[x,t],x,x]+x*D[u[x,t],x]-D[u[x,t],t]";\)\)

The function u[x,t] 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. We want to look for nonclassical symmetries. Let G=u-f[x,t]. The function G[x,t,u] is an invariant solution of the Heat equation if and only if it satisfies the invariance condition 

xse1(∂G/∂x)+xse2(∂G/∂t)+eta1(∂G/∂u)=0.

We will introduce this condition to the problem by adding it to rulesarray but with the infinitesimals xse1, xse2 and eta1 initially treated as arbitrary functions of (x,t,u) to be identified with infinitesimals later.

\(\(rulesarray = 
{"D[u[x,t],x,x]->D[u[x,t],t]-x*D[u[x,t],x]", "D[u[x,t],t]->(eta1prime[x,t,u[x,t]]/xse2prime[x,t,u[x,t]])-(xse1prime[x,t,u[x,t]]/xse2prime[x,t,u[x,t]])*D[u[x,t],x]"};\)\)

Enter the list of independent variables.

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

Enter the list of dependent variables.

dependentvariables = {"u"};

Enter the list of function and/or constant names that must be preserved when the equation is converted to generic y1[x1,x2] variables.

\(\(frozennames = {"xse1prime", "xse2prime", "eta1prime"};\)\)

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

p = 2;

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-Bäcklund 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 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.

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

equationgenericvariables

invarconditiontable

ztableofrules

Length[zdeterminingequations]

TraditionalForm[
      ColumnForm[
ToExpression[zdeterminingequations], 
        Left, Below]]

Now change eta1prime, xse1prime and xse2prime to the infinitesimals eta1, xse1 and xse2.

eta1prime[z1_, z2_, z3_] := eta1[z1, z2, z3]

xse1prime[z1_, z2_, z3_] := xse1[z1, z2, z3]

xse2prime[z1_, z2_, z3_] := xse2[z1, z2, z3]

zdeterminingequationsexpressions = 
      ToExpression[zdeterminingequations]

These are the full nonlinear determining equations. This is an extremely complex system and the only reasonable way to make progress is to look for simplifying assumptions that lead to interesting solutions. To begin set xse2=1.

zdeterminingequationssimplified = 
        zdeterminingequationsexpressions /. 
{xse2[z1, z2, 
                z3] -> 
              1, 
D[xse2[z1, z2, z3], z1] -> 
              0, 
D[xse2[z1, z2, z3], z2] -> 
              0, 
D[xse2[z1, z2, z3], z3] -> 
              0, 
D[xse2[z1, z2, z3], z1, z1] -> 
              0, 
D[xse2[z1, z2, z3], z1, z2] -> 
              0, 
D[xse2[z1, z2, z3], z1, z3] -> 
              0, 
D[xse2[z1, z2, z3], z2, z2] -> 
              0, 
D[xse2[z1, z2, z3], z2, z3] -> 
              0, 
D[xse2[z1, z2, z3], z3, z3] -> 
              0};

TraditionalForm[
      ColumnForm[
zdeterminingequationssimplified, Left, 
        Below]]

zdeterminingequationssimplified1 = 
        zdeterminingequationssimplified /. 
{D[
                xse1[z1, z2, z3], z3, z3] -> 0};

TraditionalForm[
      ColumnForm[
zdeterminingequationssimplified1, Left, 
        Below]]

This is still a fairly complex nonlinear system for xse1 and eta1. To make further progress let eta1=0 and D[xse1,z3]=0.

zdeterminingequationssimplified2 = 
        zdeterminingequationssimplified1 /. 
{eta1[z1, z2, 
                z3] -> 
              0, 
D[eta1[z1, z2, z3], z1] -> 
              0, 
D[eta1[z1, z2, z3], z2] -> 
              0, 
D[eta1[z1, z2, z3], z3] -> 
              0, 
D[eta1[z1, z2, z3], z1, z1] -> 
              0, 
D[eta1[z1, z2, z3], z1, z2] -> 
              0, 
D[eta1[z1, z2, z3], z1, z3] -> 
              0, 
D[eta1[z1, z2, z3], z2, z2] -> 
              0, 
D[eta1[z1, z2, z3], z2, z3] -> 
              0, 
D[eta1[z1, z2, z3], z3, z3] -> 
              0, 
D[xse1[z1, z2, z3], z3] -> 
              0, 
D[xse1[z1, z2, z3], z3, z1] -> 
              0, 
D[xse1[z1, z2, z3], z3, z2] -> 
              0, 
D[xse1[z1, z2, z3], z3, z3] -> 
              0};

TraditionalForm[
      ColumnForm[
zdeterminingequationssimplified2, Left, 
        Below]]

Finally with a series of simplified assumptions we have pared the determining equations to a simple result although other solutions to the full nonlinear system are certainly possible. The final PDE for xse1 can be simplified further.  Let xse1=-x+g[s], where s=x+Exp[-t]. The first equation above reduces to
g''-2gg'+sg'+g=0. Mathematica can solve this nonlinear ODE

DSolve[g''[s] - 2  g[s]*g'[s] + s*g'[s] + g[s] == 0, g[s], 
      s]

?? Erfi

The Fokker-Planck equation admits the rather complex nonclassical symmetry,

xse1=-x+g[x+Exp[-t]]
xse2=1
eta1=0

where g[s] is given above.

MaxMemoryUsed[ ]