Artifact 051fe7d3c0f86ba99a1624f63fafe7df9b2408e541221b97c98ce112b1ad9738:
- Executable file
r36/src/applysym.red
— part of check-in
[f2fda60abd]
at
2011-09-02 18:13:33
on branch master
— Some historical releases purely for archival purposes
git-svn-id: https://svn.code.sf.net/p/reduce-algebra/code/trunk/historical@1375 2bfe0521-f11c-4a00-b80e-6202646ff360 (user: arthurcnorman@users.sourceforge.net, size: 39561) [annotate] [blame] [check-ins using] [more...]
module applysym; load_package crack, crackapp; create!-package('(applysym),'(contrib crack)); fluid '(print_ fname_ nfct_ time_ facint_ logoprint_); lisp flag('(yesp),'boolean)$ lisp(logoprint_:=t)$ symbolic operator freeoflist$ %---------------------------- algebraic procedure applysym(problem,symtry); % problem ... {{equations},{functions},{variables}} % Symtry ... {{xi_..=..,..,eta_..=..,..},{constants in first list}} begin scalar genlist,con,e1,e2,h1,h2,h3,h4,h5,h6,h7,modus,u,v,xlist,ylist,n, cop1,cop2,symanz,oldsol,oldmodus,trafoprob,altlogo$ clear sy_,sym_; array sym_(length(second symtry)); symbolic put('sy_,'simpfn,'simpiden)$ symbolic put('ff,'simpfn,'simpiden)$ symbolic put('ffi,'simpfn,'simpiden)$ ylist := maklist second problem; xlist := maklist third problem; con:=second symtry; symanz:=0; for each e1 in con do % i.e. for all symmetries do: if freeoflist(e1,xlist) and freeoflist(e1,ylist) then %------ no pseudo-Lie-symm. <<genlist:=sub(e1=1,first symtry); %------ calculate symmetry for each el2 in con do if el1 neq el2 then genlist:=sub(el2=0,genlist); symanz:=symanz+1; sym_(symanz):=genlist >>; repeat << % Application oldmodus:=modus; repeat << lisp <<terpri()$ write"Do you want to find similarity and symmetry variables ", "(enter `1;')";terpri()$ write"or generalize a special solution with new parameters ", "(enter `2;')";terpri()$ write"or exit the program ", "(enter `;')";terpri() >>; modus:=aread(nil) >> until (modus=1) or (modus=2) or (modus=nil); if modus neq nil then << % Preparing a combination of symmetries if symanz=1 then genlist:=sym_(1) else << for n:=1:symanz do << write"---------------------- The ",n,". symmetry is:"$ for each e2 in sym_(n) do if rhs e2 neq 0 then write e2 >>; write"----------------------"$ repeat <<lisp <<terpri()$ write "Which single symmetry or linear combination of symmetries"$ terpri()$write"do you want to apply? "$ terpri()$ write"Enter an expression with `sy_(i)' for the i'th ", "symmetry. Terminate input with `$' or `;'."$ terpri() >>$ h1:=aread(nil); if h1 then << for each h2 in xlist do if h1 then if df(h1,h2) neq 0 then h1:=nil; for each h2 in ylist do if h1 then if df(h1,h2) neq 0 then h1:=nil; if h1=nil then lisp << terpri(); write"The coefficients of the sy_(i) must be constant, ", "i.e. numbers or constants";terpri() >> >> >> until h1; genlist:={}; cop1:=sym_(1); while cop1 neq {} do << h6:=lhs first cop1;cop1:=rest cop1; genlist:=cons(h6 = 0, genlist) >>; genlist:=reverse genlist; for h2:=1:symanz do << h3:=coeffn(h1,sy_(h2),1); if h3 neq 0 then << cop1:=genlist;cop2:=sym_(h2); genlist:={}; while cop1 neq {} do << h4:=first cop1; cop1:=rest cop1$ h5:=first cop2; cop2:=rest cop2$ h6:=lhs h4; genlist:=cons(h6 = rhs h4 + h3*(rhs h5),genlist) >>; genlist:=reverse genlist >> >> >>; write"The symmetry to be applied in the following is "; write genlist; if modus=1 then << write"Enter the name of the new dependent variable(s):"; u:=aread(nil); write"Enter the name of the new independent variable(s):"; v:=aread(nil); altlogo:=logoprint_; logoprint_:=nil; trafoprob:=similarity(problem,genlist,{},u,v); logoprint_:=altlogo >> else << if length h1 < 2 then << lisp <<terpri()$ write"What shall the name of the new constant parameter", " be? ";terpri()>>$ h2:=aread(nil) >>; repeat << lisp <<terpri()$ write"Enter the solution to be generalized in form of an ", "expression, which vanishes";terpri()$ write"or in form of an equation `... = ...' "; if oldsol=nil then write":" else << terpri()$ write"or enter semicolon `;' to work on the solution ", "specified before:" >>;terpri() >>; h3:=aread(nil); if h3 neq nil then oldsol:=h3 >> until oldsol neq nil; h3:=newparam(oldsol,genlist,h2); if h3 neq nil then oldsol:=h3 >> >> >> until modus=nil; clear sym_; return if oldmodus=1 then trafoprob else if oldmodus=2 then oldsol else nil end$ % of ApplySym %---------------------------- algebraic procedure storedepend(xlist)$ % stores the dependencies of the elements of xlist in the returned % list and clears the dependencies begin scalar q,e1,e2$ return for each e1 in xlist collect << q:=fargs e1; for each e2 in q do nodepend e1,e2; cons(e1,q)>> end$ % of storedepend %---------------------------- algebraic procedure restoredepend(prev_depend)$ % assigns the dependencies stored in the argument begin scalar q,s,x; for each s in prev_depend do <<q:=first s; s:=rest s; for each x in s do depend q,x>> end$ % of restoredepend %---------------------------- algebraic procedure newparam(oldsol,genlist,u_)$ % u_ is the name of the new parameter begin scalar h1,h2,h20,h3,h30,h4,vari,pde,e1,printold,clist,newsol, oldsol_ex,prev_depend; vari:=makepde(genlist,u_)$ pde:=first vari; vari:=rest vari; % is oldsol an invariant? oldsol_ex:=equ_to_expr(oldsol)$ % oldsol as vanishing expression prev_depend:=storedepend(vari)$ h2:=sub(u_=oldsol_ex,pde); if h2 neq 0 then << h1:=solve(oldsol_ex,vari); if h1 neq {} then << h1=first h1; for each h3 in h1 do if freeof(h3,arbcomplex) then h2:=sub(h3,h2) >> >>; restoredepend(prev_depend)$ if 0=h2 then return lisp <<write"The special solution to be generalized is an invariant ", "with respect to";terpri()$ write"this symmetry, therefore no generalization is possible."; terpri()$ for each h1 in fargs u_ do nodepend u_,h1; algebraic oldsol >>; pde:=pde-1; h1:= quasilinpde1(pde,u_,vari)$ if h1 neq {} then << h1:=first h1; % h2 is expressing the constants in terms of u_ and the xlist,ylist clist:={}; lisp<<printold:=print_; print_:=nil>>; h2:= for each e1 in h1 collect <<h3:=lisp(newfct(fname_,nil,nfct_))$ lisp(nfct_:=add1 nfct_)$ clist:=cons(h3,clist); h3 = e1>>; lisp(print_:=printold); h20:=sub(u_=0,h2); h3:=solve(h2,vari); if h3 neq {} then << h3:=first h3; h30:=sub(h20,h3); write"The substitutions to generalize the solution are: "$ for each h4 in h30 do write h4$ newsol:=sub(h30,oldsol_ex); % newsol:=second dropredundant(newsol, % alle Konstanten und % Funktionen die nicht in vari sind, % vari)$ lisp <<write"The new solution"; % if length algebraic h >2 then write"s are:" else write" is:"; % terpri()$ >>; write"0 = ",newsol >> >>; for each h1 in fargs u_ do nodepend u_,h1; return newsol end$ % of NewParam %---------------------------- %symbolic operator aread$ %symbolic procedure aread(prompt)$ %begin scalar val, !*echo; % Don't re-echo tty input % rds nil; wrs nil$ % Switch I/O to terminal % terpri()$ % if null prompt then prompt:= "Input:"; prin2 prompt; % val := xread(nil)$ % if ifl!* then rds cadr ifl!*$ % Resets I/O streams % if ofl!* then wrs cdr ofl!*$ % return val %end$ %---------------------------- symbolic operator einfachst$ symbolic procedure einfachst(a,x)$ % a is an algebraic list and the element where x appears, but % appears simplest, is found begin scalar el1,el2,hp,printold; printold:=print_; print_:=nil; hp:=10000; a:=cdr a; while a do << el2:=car a;a:=cdr a; if not freeof(el2,x) and ((not el1 ) or (el2 = x ) or <<if not polyp(el2,cons(x,nil)) then nil else <<coeff1(el2,x,nil)$ if hipow!*<hp then <<hp:=hipow!*; t>> else nil >> >> ) then el1:=el2 >>; print_:=printold; return el1 end$ %---------------------------- algebraic procedure transdf(y,yslist,vlist,indxlist)$ begin scalar m,n,e1,dfy; return if indxlist={} then sub(yslist,y) else << m:=first indxlist; n:=0; dfy:=transdf(y,yslist,vlist,rest indxlist); for each e1 in vlist sum << n:=n+1; df(dfy,e1)*dv!/dx(n,m) >> >> end$ % of TransDf %---------------------------- algebraic procedure transderiv(yik,yslist,vlist)$ begin scalar indxlist,y,l1,l2; indxlist:=lisp cons('list,combidif(yik))$ return transdf(first indxlist,yslist,vlist,rest indxlist) end$ % of TransDeriv %---------------------------- algebraic procedure detrafo(eqlist,yslist,xslist,ulist,vlist)$ % Transformations of all orders are performed (point-, contact-,...) % but only x-derivatives of y's are transformed. To transform other % any other derivatives, subdiff1 must be extended to include all other % occuring x-derivatives. begin scalar avar,nvar,detpd,n,m,ordr,e1,e2,e3,sb; m:=length(xslist); n:=length(yslist)+m; clear dyx!/duv,dv!/dx; matrix dyx!/duv(n,n); matrix dv!/dx(m,m); avar:=append(yslist,xslist); nvar:=append(ulist,vlist); n:=0; for each e1 in avar do << n:=n+1;m:=0; for each e2 in nvar do << m:=m+1; dyx!/duv(m,n):=df(rhs e1,e2) >> >>; detpd:=det(dyx!/duv); %write"detpd=",detpd; if detpd=0 then return <<write"The proposed transformation is not regular!";{}>>; clear dyx!/duv; ordr:=0; for each e1 in eqlist do for each e2 in yslist do <<n:=totdeg(e1,lhs e2); if n>ordr then ordr:=n>>; sb:=subdif1(for each e1 in xslist collect lhs e1, for each e1 in yslist collect lhs e1,ordr); % computation of Dv/Dx:=(Dx/Dv)^(-1) n:=0; for each e1 in xslist do << n:=n+1;m:=0; for each e2 in vlist do << m:=m+1; dv!/dx(n,m):=df(rhs e1,e2) >> >>; dv!/dx:=dv!/dx**(-1); % Substitution of all derivatives for each e1 in sb do for each e2 in e1 do << if not freeof(eqlist,lhs e2) then << % which function is to be differentiated wrt. which variable eqlist:=sub(lhs e2=transderiv(rhs e2,yslist,vlist),eqlist) >> >>; clear dv!/dx; return sub(xslist,sub(yslist,eqlist))$ end$ % of DeTrafo %---------------------------- algebraic procedure grouping(el1,el2,xlist,ylist,nx,ny)$ begin scalar h,el3,xslist,yslist$ %------- Grouping the new variables to ulist and vlist h:={}; xslist:={}; % list of expressions to calculate new indep. var. yslist:={}; % list of expressions to calculate new dep. var. %---- at first the obvious allocations for each el3 in el1 do %-- all similarity variables if freeoflist(el3,ylist) then xslist:=cons(el3,xslist) else if freeoflist(el3,xlist) then yslist:=cons(el3,yslist) else h:=cons(el3,h); %---- now the symmetry variable if freeoflist(el2,ylist) or (length(yslist) = ny) then xslist:=cons(el2,xslist) else if freeoflist(el2,xlist) or (length(xslist) = nx) then yslist:=cons(el2,yslist) else xslist:=cons(el2,xslist); %---- now the remaining cases for each el3 in h do if length(yslist) < ny then yslist:=cons(el3,yslist) else xslist:=cons(el3,xslist); return {xslist,yslist} end$ % of grouping %---------------------------- algebraic procedure rename_u_(xslist,yslist,el2,u_,u,v)$ begin scalar i,vlist,ulist,el3,h,smv$ %---- Renaming the u_ to ui in yslist and to vi in xslist i:=0; vlist:={}; xslist:=for each el3 in xslist collect <<i:=i+1; if length xslist>1 then h:=mkid(v,i) else h:=v; vlist:=cons(h,vlist); if el3=el2 then smv:=h; sub(u_=h,el3) >>; i:=0; ulist:={}; yslist:=for each el3 in yslist collect <<i:=i+1; if length yslist>1 then h:=mkid(u,i) else h:=u; ulist:=cons(h,ulist); if el3=el2 then smv:=h; sub(u_=h,el3) >>; return {xslist,yslist,reverse vlist,reverse ulist,smv} end$ % of rename_u_ %---------------------------- algebraic procedure solve_for_old_var(xslist,yslist,xlist,ylist,nx,ny)$ begin scalar h1,h2$ %---- Solve for old variables h1:=nil; h2:=solve(append(yslist,xslist),append(xlist,ylist)); if h2={} then h1:=t else h2:=first h2; %--- possibly other solutions if list neq lisp(car algebraic h2) then h1:=t else if length(h2)<(nx+ny) then el2:=t; if h1 then repeat lisp <<write"The algebraic system ",append(xslist,yslist), " could not be solved for ",append(xlist,ylist),"."; write"Please enter the solution in form of a list {", reval algebraic first xlist,"=...,...", reval algebraic first ylist,"=...,...} or enter a ", "semicolon ; to end this investigation:"; algebraic(h2:=aread(nil)) >> until h2=nil or ( lisp(pairp algebraic h2) and (list=lisp(car algebraic h2)) and (length(h2)=(nx+ny)) ) else <<lisp<<terpri()$ write"The suggested solution of the algebraic system which will"; terpri()$ write"do the transformation is: "; terpri() >>; write h2; if yesp "Is the solution ok?" then else lisp << write"Please enter the solution in form of a list {", reval algebraic first xlist,"=...,...", reval algebraic first ylist,"=...,...} or enter a ", "semicolon ; to end this investigation:"; algebraic(h2:=aread(nil)) >> >>; return h2 end$ % of solve_for_old_var$ %---------------------------- algebraic procedure switch_r_s(h2,smv,ylist,u,v)$ begin scalar xslist,yslist,el3,h$ %---- Exchange of dependent and independent variables xslist:={}; yslist:={}; for each el3 in h2 do if freeof(ylist,lhs el3) then xslist:=cons(el3,xslist) else yslist:=cons(el3,yslist); lisp <<terpri()$ write"In the intended transformation shown above", " the dependent ";terpri()$ if length yslist>2 then write"variables are the ",reval algebraic u,"i and " else write"variable is ",reval algebraic u," and "; if length xslist>2 then write"the independent variables are the ", reval algebraic v,"i." else write"the independent variable is ",reval algebraic v,"."; terpri()$ write "The symmetry variable is ",reval algebraic smv,", i.e. the ", "transformed expression";terpri(); write"will be free of ",reval algebraic smv,"."; >>; h:=if yesp "Is this selection of dependent and independent variables ok?" then nil else << lisp <<write"Please enter a list of substitutions. For example, to"; terpri()$ write"make the variable, which is so far call u1, to an"; terpri()$ write"independent variable v2 and the variable, which is "; terpri()$ write"so far called v2, to an dependent variable u1, "; terpri()$ write"enter: `{u1=v2, v2=u1};'">>; aread(nil) >>; if h and (h neq {}) then <<xslist:=sub(h,xslist); yslist:=sub(h,yslist); smv :=sub(h,smv)>>; return {xslist,yslist,smv} end$ % of switch_r_s %---------------------------- algebraic procedure makepde(genlist,u_)$ begin scalar h,el2,el3,vari,bv; vari:={}; return cons( num for each el2 in genlist sum <<h:=lhs el2; h:=lisp << el3:=explode reval algebraic h; bv:=t; while bv do << if car el3 ='!_ then bv:=nil; el3:=cdr el3 >>; intern compress el3 >>; depend u_,h; vari:=cons(h,vari); (rhs el2) * df(u_,h) >>, vari) end$ % of makepde %---------------------------- algebraic procedure totdeglist(eqlist,ylist)$ begin scalar n,ordr,e1,e2; ordr:=0; for each e1 in eqlist do for each e2 in ylist do <<n:=totdeg(e1,e2); if n>ordr then ordr:=n>>; return ordr end$ % of totdeglist %---------------------------- algebraic procedure similarity(problem,genlist,con,u,v)$ % con ... the free constants/functions in the general symmetry % u ... the name of the new independent variables % v ... the name of the new dependent variables begin scalar vari,pde,el1,el2,el3,el4,copgen,symvarfound, trans1,trans2,i,j,h,h2,n,denew,xlist,ylist, eqlist,ulist,vlist,nx,ny,xslist,yslist,smv, trafoprob; cpu:=lisp time()$ gc:=lisp gctime()$ %--------- extracting input data eqlist:=maklist first problem; ylist :=maklist second problem; ny:=length ylist; xlist :=maklist third problem; nx:=length xlist; trafoprob:={problem,nil}; % to be returned if trafo not possible problem:=nil; eqlist:=for each el1 in eqlist collect equ_to_expr el1; if length eqlist > 1 then eqlist:=desort eqlist; %--------- initial printout lisp( % if print_ and logoprint_ then << % write "-----------------------------------------------", % "---------------------------"$ terpri()$terpri()$ % write "This is SIMILARITY - a program for calculating similarity"; % write % "variables and performing corresponding point-transformations"; % write "for given infinitesimal symmetries."; % terpri() % >> else if print_ then terpri()); %--------- initializations ordr:=totdeglist(eqlist,ylist)$ vari:=append(ylist,xlist); for each el1 in xlist do for each el2 in ylist do nodepend el2,el1; lisp( if print_ then << write "The ODE/PDE (-system) under investigation is :";terpri()$ for each el1 in cdr eqlist do algebraic write"0 = ",el1; terpri()$write "for the function(s) : "; fctprint( cdr reval algebraic ylist);write"."; terpri()$terpri() >>); lisp( if print_ then << if length ylist >2 then % not >1 because of alg. list in symb. mode write"It will be looked for new dependent variables ",u,"i " else write"It will be looked for a new dependent variable ",u; terpri()$ if length xlist >2 then write"and independent variables ",v,"i" else write"and an independent variable ",v; write" such that the transformed"; terpri()$ write"de(-system) does not depend on ",u; if length ylist >2 then write"1"; write" or ",v$ if length xlist >2 then write"1"; write"."; terpri() >>); % for each el1 in con do % i.e. for all symmetries do: % if freeoflist(el1,xlist) and % freeoflist(el1,ylist) then %------ no pseudo-Lie-symm. % <<copgen:=sub(el1=1,genlist); %------ calculate symmetry % for each el2 in con do % if el1 neq el2 then copgen:=sub(el2=0,copgen); copgen:=genlist; % write"The symmetry now under investigation is:"; % for each el1 in copgen do write el1; %---------- formulate the PDE for the similarity variables pde:=first makepde(copgen,u_); %--------- find similarity variable trans2 :={}; lisp<<terpri()$write"1. Determination of the similarity variable"; if nx+ny>2 then write"s">>; trans1 := quasilinpde1(pde,u_,vari); % for the similarity variable if trans1 neq {} then << %-------------- Determining the similarity variables ui_ i:=0; trans1:=for each el1 in trans1 collect <<i:=i+1; h:=length(genlist)-1; if h=1 then % one single ODE << % write"In the following 1 similarity variable U_ has to be"; % write"determined through 0 = ff where ff is an"$ % write"arbitrary function of arguments given in the ", % "following list:"$ % write el1; el2:=num(first el1 - second el1); if freeof(el2,u_) then el2:=num(first el1 - 2*second el1); lisp<<write"A suggestion for this function ff provides:"$ terpri()>>$ write"0 = ", el2$ if yesp "Do you like this choice?" then {el2} else << repeat << lisp << write"Put in an alternative expression which "$terpri()$ write"- is functionally dependent only on elements of", " ff given above and"$ terpri()$ write"- depends on U_ and if set to zero determines U_"$ terpri()>>$ h:=aread(nil)$ >> until not freeof(h,u_)$ {h} >> >> else % a PDE or a system of DEs lisp << % terpri()$ % write"Now the similarity variables U_i, i=1,...",h, % " have to be"$terpri()$ % write"determined through conditions 0 = ffi, i=1,...",h, % ", where ffi are"$terpri()$ % write"arbitrary functions ";terpri()$ % algebraic(write "ffi = ", % lisp( cons('ffi,cdr reval algebraic el1))); % terpri()$ % write"such that the functional determinant of these ", % "expressions"$ terpri()$ % write"including u_ from above taken w.r.t. "; % for each el3 in cdr algebraic ylist do write reval el3,","; % write reval cadr algebraic xlist; % for each el3 in cddr algebraic xlist do write ",",reval el3; % terpri()$ % write"must not vanish."$terpri()$ algebraic << el2:=einfachst(el1,u_); h2:={}; for each el3 in el1 do if el3 neq el2 then h2:=cons(num(el2-el3),h2) >>; write "A suggestion for these functions ffi in form of a list ", "{ff1,ff2,... } is: "$terpri()$ deprint(cdr reval algebraic h2); if yesp "Do you like this choice?" then algebraic h2 else << write"Put in an alternative list of expression which"$ terpri()$ write"- are functionally dependent only on the above ", "arguments and"$terpri()$ write"- which if set to zero determine U_i, i.e."$terpri()$ write"- the functional determinant of these expressions"$ terpri()$ write" including U_ from above taken w.r.t. ", cdr reval algebraic append(ylist,xlist)$ terpri()$ write" must not vanish."$terpri()$ algebraic(h2):=aread(nil) >> >> >>$ %--------- find symmetry variable pde:=pde-1; lisp<<terpri()$write"2. Determination of the symmetry variable">>; trans2 := quasilinpde1(pde,u_,vari); % for the symmetry variable if (length xlist=1) and (trans2={}) then for each e1 in fargs u_ do nodepend u_,e1 else <<% If no symmetry variable is found (trans2={}) then proceed % only if special solutions of PDEs are to be found with the % solution being only a function of the similarity variables if trans2={} then << % take any variable h:=reverse vari; % reverse to have not a function as symvar while 1+sub(u_=first h,pde)=0 do h:=rest h; % no similarity var. for each e1 in fargs u_ do nodepend u_,e1; h:=first h; trans2:={u_ - h}; lisp<<write"Because the correct symmetry variable was not ", "found, the program will";terpri()$ write"take ",reval algebraic h, " instead with the consequence ", "that not the whole transformed ";terpri()$ write"PDE will be free of ", reval algebraic h," but only those ", "terms without ", reval algebraic h,"-derivative";terpri()$ write"which is still of use for finding special ", reval algebraic h,"-independent solutions "; terpri()$ write"of the PDE." >> >> else << %--------------- Determining an optimal symmetry variable for each e1 in fargs u_ do nodepend u_,e1; symvarfound:=t; i:=0; trans2:=for each el1 in trans2 collect <<i:=i+1; lisp << % terpri()$ % write"In the following the symmetry variable U_", % " has to be"$terpri()$ % write"determined through a condition 0 = ff", % ", where ff is"$terpri()$ % write"an arbitrary function ";terpri()$ % algebraic(write"ff = ", % lisp( cons('ff,cdr reval algebraic el1)));terpri()$ write "A suggestion for this function ff(..) yields:";terpri() >>; h:=einfachst(el1,u_); if lisp<<h2:=reval algebraic(num h); (not pairp h2) or (car h2 neq 'plus)>> then h:=num(h+1); % if h= first el1 then % if freeof(num(h+second el1),u_) then h:=num(h+2*second el1) % else h:=num(h+ second el1) % else % if freeof(num(h+ first el1),u_) then h:=num(h+2* first el1) % else h:=num(h+ first el1); write"0 = ",h$ if yesp "Do you like this choice?" then h else << repeat << lisp << write "Put in an alternative expression which "$terpri()$ write "- is functionally dependent only on arguments of", " ff given above and"$terpri()$ write "- depends on u_ and if set to zero determines u_"$ terpri()>>$ h:=aread(nil)$ >> until not freeof(h,u_)$ h >> >> >>$ %for each el1 in trans1 do %for each el2 in trans2 do el1:=first trans1; el2:=first trans2; % << %------- Grouping the new variables to ulist and vlist yslist:=grouping(el1,el2,xlist,ylist,nx,ny)$ xslist:= first yslist; yslist:=second yslist; %---- Renaming the u_ to ui in yslist and to vi in xslist smv:=rename_u_(xslist,yslist,el2,u_,u,v)$ xslist:= first smv; yslist:=second smv; vlist := third smv; ulist:=part(smv,4); smv:=part(smv,5); % the symmetry variable %---- Solve for old variables h2:=solve_for_old_var(xslist,yslist,xlist,ylist,nx,ny); if h2 neq nil then << %---- Exchange of dependent and independent variables smv:=switch_r_s(h2,smv,ylist,u,v)$ xslist:= first smv; yslist:=second smv; smv :=third smv; %---- Doing the point transformation for each el3 in ulist do <<for each el4 in fargs el3 do nodepend el3,el4; for each el4 in vlist do %if el4 neq smv then % if new DEs without symm. var. smv depend el3,el4>>; for each el3 in ylist do for each el4 in xlist do depend el3,el4; eqlist:=detrafo(eqlist,yslist,xslist,ulist,vlist); lisp( if print_ then << terpri()$ write"The transformed equation"; if length(algebraic eqlist)>2 then write"s"; if symvarfound then write" which should be free of ",reval algebraic smv,":" else write" in which the terms without ",reval algebraic smv, "-derivative are free of ",reval algebraic smv,":"; terpri() >>); eqlist:=for each el3 in eqlist collect << el3:=factorize num el3; for each el4 in el3 product if 0=totdeglist({el4},ulist) then 1 else el4 >>; lisp deprint(cdr reval algebraic eqlist); if (length(vlist)>1) and (not freeof(vlist,smv)) then << vlist:=cons(smv,lisp(delete(reval algebraic smv, reval algebraic vlist))); if yesp "Shall the dependence on the symmetry variable be dropped?" then <<for each el3 in ulist do nodepend el3,smv; vlist:=rest vlist>>; eqlist:=for each el3 in eqlist collect << el3:=factorize num el3; for each el4 in el3 product if 0=totdeglist({el4},ulist) then 1 else el4 >> >>; trafoprob:={{eqlist,ulist,vlist},append(xslist,yslist)} >> % >> >> >>; % >>; for each el1 in xlist do for each el2 in ylist do depend el2,el1; clear ff,ffi; return trafoprob end$ % of similarity %---------------------------- symbolic procedure simplifiziere(q)$ begin scalar n,one_argument_functions_; one_argument_functions_:= '(abs acos acosd acosh acot acotd acoth acsc acscd acsch asec asecd asech asin asind asinh atan atand atanh cbrt cos cosd cosh cot cotd coth csc cscd csch exp hypot ln log log10 sec secd sech sin sind sinh sqrt tan tand tanh minus)$ return if not pairp q then q else if member(car q,one_argument_functions_) or ((car q = 'expt) and (numberp caddr q)) then simplifiziere(cadr q) else << n:=ncontent(q); if n=1 then q else simplifiziere(reval list('quotient, q, reval n)) >> end$ %---------------------------- algebraic procedure select_indep_var(clist,xlist)$ begin scalar s,scal,notfound,cq,x,xx,xanz,sanz,xcop,ok; % Find the simplest non-zero element of clist notfound:=t; xcop:=xlist; while notfound and (xcop neq {}) do << x :=first xcop ; xcop :=rest xcop ; cq:=first clist; clist:=rest clist; if cq neq 0 then << xanz:=0; for each xx in xlist do if %df(cq,xx) not freeof(cq,xx) then xanz:=xanz+1; ok:=nil; if not s then ok:=t else % to guarantee s neq nil if xanz<sanz then ok:=t else % as few dependencies as poss. if xanz=sanz then if length(cq)=1 then ok:=t else if length(scal)>1 then if part(scal,0)=plus then % if possible not a sum if (part(cq,0) neq plus) or (length(cq)<length(scal)) then ok:=t; if ok then <<s:=x; scal:=cq; sanz:=xanz>>; if scal=1 then notfound:=nil >> >>; return {s,scal} end$ % of select_indep_var %---------------------------- procedure charsyspsfi(xlist,cqlist,u,ode,denf); begin scalar h,s; h:=cqlist; cqlist:={}; while h neq {} do <<cqlist:=cons(first h*denf,cqlist);h:=rest h>>; cqlist:=cons(-ode*denf,cqlist); xlist:=cons(u,reverse xlist); s:=lisp gensym(); for each h in xlist do depend h,s; h:=psfi(cqlist,xlist,s); write "soln=",h; for each h in xlist do nodepend h,s; return h end$ % of charsyspsfi %---------------------------- algebraic procedure charsyscrack(xlist,cqlist,s,scal,u,ode)$ % formulation and solution of the characteristic ODE-system begin scalar lcopy,x,cs,flist,soln,printold,timeold,facintold, e1,e2,e3,e4,n,nmax,dff; % formulating the ode/pde - system . lcopy := cqlist ; cs := {} ; for each x in xlist do << cq := first lcopy ; lcopy := rest lcopy ; if x neq s then << depend x,s; cs := .(scal*df(x,s)-cq,cs) >> >>; if s neq u then <<flist:={u}; depend u,s; cs := .(scal*df(u,s)+ode,cs) >> else flist:={}; for each x in xlist do if s neq x then flist:=.(x,flist); % solving the odes using crack. if lisp(print_ neq nil) then lisp << write "The equivalent characteristic system:";terpri(); deprint cdr algebraic cs$ %terpri()$ write "for the functions: "; fctprint( cdr reval algebraic flist);write"."; printold:=print_; print_:=nil; timeold := time_; time_ :=nil; facintold:=facint_;facint_:=t; >>; % lisp(print_:=100); %############ soln:=crack(cs,flist,flist,{}); lcopy:={}; for each x in soln do << e1:=first x; if (e1 neq {}) and (length e1 + length second x = length flist) then % all remaining conditions are algebraic (non-differential) in all % the functions of flist? <<e2:={}; e3:={}; for each e4 in third x do if freeof(flist,e4) then e3:=e4 . e3 else e2:=e4 . e2; if (length cs) = (length e3) then << % sufficiently many integrations done for each e4 in e2 do for each e5 in e1 do if lisp(not freeof(lderiv(reval algebraic e5, reval algebraic e4, list(reval algebraic s)),'df)) then <<e2:={};e1:={}>>; if e2 neq {} then << % It may be possible that derivatives of unsolved functions % occur in the expressions of the evaluated functions: second x nmax:=0; for each e4 in e2 do << % for each unsolved function for each e5 in second x do << % for each solved expression n:=lisp(lderiv(reval algebraic rhs e5,reval algebraic e4, list(reval algebraic s))); if (n neq 0) and (n neq nil) then if lisp length algebraic n = 3 then n:=1 else n:=lisp cadddr algebraic n; if n>nmax then nmax:=n; >>; if nmax>0 then << % adding further conditions e5:=e1; while freeof(first e5,e4) do e5:=rest e5; e5:=first e5; dff:=e4; for n:=1:nmax do << e5 :=df(e5,s); e1:= e5 . e1; dff:=df(dff,s); e3:=dff . e3 >> >> >>; lcopy:=cons({append(second x,e1),e3},lcopy); >> >> >> else if (first x = {}) and (length cs = length third x) then lcopy:=cons({second x,third x},lcopy) >>; lisp <<print_:=printold; time_:=timeold; facint_:=facintold>>; return if lcopy={} then <<for each x in flist do nodepend x,s; {}>> else s . lcopy % { {{x=..,y=..,u=..,..,0=..},{df(z,s),df(z,s,2),..,c1,c2,c3,..}},..} % df(z,s,..) only if df(z,s,..) turns up in x, y, u. .. . end$ % of charsyscrack %---------------------------- algebraic procedure quasilinpde1(pde,u_,vari)$ begin scalar trans1,e1,e2,q; trans1 := quasilinpde(pde,u_,vari); % for the similarity variable if trans1={} then << write"The program was not able to find the general solution ", "of the PDE: ",pde," for the function ",u_,"."; lisp << write"Please enter either only a semicolon if no solution ", "is known or enter "$terpri()$ write"the solution of the PDE in form ", "of a list {A1,A2,...} where ";terpri()$ write"the Ai are algebraic expressions in ", cdr reval algebraic cons(u_,vari);terpri()$ write"such that any function ff(A1,A2,...) which is not ", "independent of `",u_;write"'";terpri()$ write"determines a solution `",u_,"' of the PDE through 0=ff: " >>; trans1:={aread(nil)}; if trans1={nil} then trans1:={} >>; return trans1 end$ % of quasilinpde1 %---------------------------- algebraic procedure quasilinpde(f,u,xlist); % f ... PDE, u ... unknown function, xlist ... variable list begin scalar i, q, qlist, cq, cqlist, ode, soln, tran, xcop, s, s1, x, xx, scal, qlin, prev_depend; symbolic put('ff,'simpfn,'simpiden)$ % if lisp(print_ neq nil) then write"The quasilinear PDE: 0 = ",f,"."; % changing the given pde into a quasi-linear ode . i := 0 ; ode := f; qlist := {}; cqlist:={}; for each x in xlist do <<i := i+1 ; q := mkid(p_,i) ; qlist := .(q,qlist) ; f := sub(df(u,x) = q , f) ; cq:=df(f,q); cqlist:=.(df(f,q),cqlist) ; ode := ode - df(u,x)*df(f,q) ; nodepend u, x ; >> ; qlist := reverse qlist ; cqlist := reverse cqlist ; % checking for linearity. qlin:=t; for each cq in cqlist do for each q in qlist do if df(cq,q) neq 0 then qlin:=nil; if not qlin then return {}; % soln:=charsyspsfi(xlist,cqlist,u,ode,den f); % Determination of the independent variable for the ODEs pcopy:=cons(-ode,cqlist);xcop:=cons(u,xlist); scal:=select_indep_var(pcopy,xcop)$ s1:=first scal; prev_depend:=storedepend(xlist)$ soln:=charsyscrack(xlist,cqlist,s1,second scal,u,ode); if soln={} then << % try all other coefficients as factors repeat << repeat <<s :=first xcop ;xcop :=rest xcop; scal:=first pcopy;pcopy:=rest pcopy>> until (pcopy={}) or ((scal neq 0) and (s neq s1)); if (s neq s1) and (scal neq 0) then << lisp <<write"New attempt with a different independent variable"; terpri()>>$ soln:=charsyscrack(xlist,cqlist,s,scal,u,ode) >> >> until (soln neq {}) or (xcop={}) >>; % solving for the constants(eg..c1,c2 etc.) and put them in a % linear form. % The first element of soln is the ODE-variable tran:={}; if soln neq {} then << s1:=first soln; for each s in rest soln do <<cq:=first solve(first s,second s); tran:=.(<<x:={}; for each xx in cq do if lisp atom algebraic lhs xx then x:=.(lisp <<q:=reval algebraic rhs xx; simplifiziere q >> ,x); x >>, % intsimplify(q,cons(reval algebraic u, % cdr reval algebraic xlist)) >>, tran); >>; for each s in xlist do if s neq s1 then nodepend s,s1 >>; for each x in xlist do depend u,x; if lisp(print_ neq nil) then if tran neq {} then << write"The general solution of the PDE is given through"; for each x in tran do write"0 = ", lisp( cons('ff,cdr reval algebraic x)); if length tran>1 then write"with arbitrary function(s) ff(..)." else write"with arbitrary function ff(..)." >>; % restoring the dependencies of the variables of the PDE restoredepend(prev_depend)$ return tran; end$ % of quasilinpde endmodule; end;