applysym.red 28 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280281282283284285286287288289290291292293294295296297298299300301302303304305306307308309310311312313314315316317318319320321322323324325326327328329330331332333334335336337338339340341342343344345346347348349350351352353354355356357358359360361362363364365366367368369370371372373374375376377378379380381382383384385386387388389390391392393394395396397398399400401402403404405406407408409410411412413414415416417418419420421422423424425426427428429430431432433434435436437438439440441442443444445446447448449450451452453454455456457458459460461462463464465466467468469470471472473474475476477478479480481482483484485486487488489490491492493494495496497498499500501502503504505506507508509510511512513514515516517518519520521522523524525526527528529530531532533534535536537538539540541542543544545546547548549550551552553554555556557558559560561562563564565566567568569570571572573574575576577578579580581582583584585586587588589590591592593594595596597598599600601602603604605606607608609610611612613614615616617618619620621622623624625626627628629630631632633634635636637638639640641642643644645646647648649650651652653654655656657658659660661662663664665666667668669670671672673674675676677678679680681682683684685686687688689690691692693694695696697698699700701702703704705706707708709710711712713714715716717718719720721722723724725726727728729730731732733734735736737738739740741742743744745746747748749750751752753754755756757758759760761762763764765766767768769770771772773774775776777778779780781782783784785786787788789790791792793794795796797798799800801802803804805806807808809810811812813814815816817818819820821822823824825826827828829830831832833834835836837838839840841842843844845846847848849850851852853854855856857858859860861862863864865866867868869870871872873874875876877878879880881882883884
  1. %********************************************************************
  2. % *
  3. % The program APPLYSYM for applying point-symmetries which are, *
  4. % for example computed by the program LIEPDE. It also can be used *
  5. % for solving quasilinear first order PDEs using QUASILINPDE. *
  6. % *
  7. % Author: Thomas Wolf *
  8. % Date: summer 1995 *
  9. % *
  10. %********************************************************************
  11. symbolic fluid '(print_ logoprint_ nfct_ fname_ time_ facint_
  12. adjust_fnc safeint_ freeint_ odesolve_)$
  13. lisp flag('(yesp),'boolean)$
  14. lisp(logoprint_:=t)$
  15. symbolic operator freeoflist$
  16. symbolic operator termxread$
  17. %----------------------------
  18. symbolic fluid '(tr_as)$
  19. lisp(tr_as:=t)$
  20. algebraic procedure ApplySym(problem,Symtry);
  21. % problem ... {{equations},{functions},{variables}}
  22. % Symtry ... {{xi_..=..,..,eta_..=..,..},{constants in first list}}
  23. begin
  24. scalar genlist,con,e1,e2,h1,h2,h3,h4,h5,h6,h7,modus,u,v,xlist,ylist,n,
  25. cop1,cop2,symanz,oldsol,oldmodus,trafoprob,altlogo$
  26. backup_reduce_flags()$
  27. clear sy_,sym_;
  28. array sym_(length(second Symtry));
  29. symbolic put('sy_,'simpfn,'simpiden)$
  30. symbolic put('ff,'simpfn,'simpiden)$
  31. symbolic put('ffi,'simpfn,'simpiden)$
  32. ylist := maklist second problem;
  33. xlist := maklist third problem;
  34. con:=second Symtry;
  35. symanz:=0;
  36. for each e1 in con do % i.e. for all symmetries do:
  37. if freeoflist(e1,xlist) and
  38. freeoflist(e1,ylist) then %------ no pseudo-Lie-symm.
  39. <<genlist:=sub(e1=1,first Symtry); %------ calculate symmetry
  40. for each el2 in con do
  41. if el1 neq el2 then genlist:=sub(el2=0,genlist);
  42. symanz:=symanz+1;
  43. sym_(symanz):=genlist
  44. >>;
  45. repeat <<
  46. % Application
  47. oldmodus:=modus;
  48. repeat <<
  49. lisp <<terpri()$
  50. write"Do you want to find similarity and symmetry variables ",
  51. "(enter `1;')";terpri()$
  52. write"or generalize a special solution with new parameters ",
  53. "(enter `2;')";terpri()$
  54. write"or exit the program ",
  55. "(enter `;')";terpri()
  56. >>;
  57. modus:=termxread()
  58. >> until (modus=1) or (modus=2) or (modus=nil);
  59. if modus neq nil then <<
  60. % Preparing a combination of symmetries
  61. if symanz=1 then genlist:=sym_(1)
  62. else <<
  63. for n:=1:symanz do <<
  64. write"---------------------- The ",n,". symmetry is:"$
  65. for each e2 in sym_(n) do
  66. if rhs e2 neq 0 then write e2
  67. >>;
  68. write"----------------------"$
  69. repeat
  70. <<lisp
  71. <<terpri()$
  72. write"Which single symmetry or linear combination of symmetries"$
  73. terpri()$write"do you want to apply? "$
  74. terpri()$write"Enter an expression with `sy_(i)' for the i'th ",
  75. "symmetry. Terminate input with `$' or `;'."$
  76. terpri()
  77. >>$
  78. h1:=lisp(reval termxread());
  79. if h1 then <<
  80. for each h2 in xlist do if h1 then if df(h1,h2) neq 0 then h1:=nil;
  81. for each h2 in ylist do if h1 then if df(h1,h2) neq 0 then h1:=nil;
  82. if h1=nil then lisp <<
  83. terpri();write"The coefficients of the sy_(i) must be constant, ",
  84. "i.e. numbers or constants";terpri()
  85. >>
  86. >>
  87. >> until h1;
  88. genlist:={};
  89. cop1:=sym_(1);
  90. while cop1 neq {} do <<
  91. h6:=lhs first cop1;cop1:=rest cop1;
  92. genlist:=cons(h6 = 0, genlist)
  93. >>;
  94. genlist:=reverse genlist;
  95. for h2:=1:symanz do <<
  96. h3:=coeffn(h1,sy_(h2),1);
  97. if h3 neq 0 then <<
  98. cop1:=genlist;cop2:=sym_(h2);
  99. genlist:={};
  100. while cop1 neq {} do <<
  101. h4:=first cop1; cop1:=rest cop1$
  102. h5:=first cop2; cop2:=rest cop2$
  103. h6:=lhs h4;
  104. genlist:=cons(h6 = rhs h4 + h3*(rhs h5),genlist)
  105. >>;
  106. genlist:=reverse genlist
  107. >>
  108. >>
  109. >>;
  110. write"The symmetry to be applied in the following is ";
  111. write genlist;
  112. write"Terminate the following input with `$' or `;'."$
  113. if modus=1 then <<
  114. write"Enter the name of the new dependent variable(s):";
  115. u:=termxread();
  116. write"Enter the name of the new independent variable(s):";
  117. v:=termxread();
  118. altlogo:=logoprint_;
  119. logoprint_:=nil;
  120. trafoprob:=similarity(problem,genlist,{},u,v);
  121. logoprint_:=altlogo
  122. >> else <<
  123. if length h1 < 2 then <<
  124. lisp <<terpri()$
  125. write"What shall the name of the new constant parameter",
  126. " be? ";terpri()>>$
  127. h2:=termxread()
  128. >>;
  129. repeat <<
  130. lisp <<terpri()$
  131. write"Enter the solution to be generalized in form of an ",
  132. "expression, which vanishes";terpri()$
  133. write"or in form of an equation `... = ...' ";
  134. if oldsol=nil then write":" else <<
  135. terpri()$
  136. write"or enter semicolon `;' to work on the solution ",
  137. "specified before:"
  138. >>;terpri()
  139. >>;
  140. h3:=termxread();
  141. if h3 neq nil then oldsol:=h3
  142. >> until oldsol neq nil;
  143. h3:=NewParam(oldsol,genlist,h2);
  144. if h3 neq nil then oldsol:=h3
  145. >>
  146. >>
  147. >> until modus=nil;
  148. clear sym_;
  149. recover_reduce_flags()$
  150. return if oldmodus=1 then trafoprob else
  151. if oldmodus=2 then oldsol else nil
  152. end$ % of ApplySym
  153. %----------------------------
  154. algebraic procedure NewParam(oldsol,genlist,u_)$
  155. % u_ is the name of the new parameter
  156. begin
  157. scalar h1,h2,h20,h3,h30,h4,vari,pde,e1,printold,clist,newsol,
  158. oldsol_ex,prev_depend;
  159. vari:=makepde(genlist,u_)$
  160. pde:=first vari;
  161. vari:=rest vari;
  162. % is oldsol an invariant?
  163. oldsol_ex:=equ_to_expr(oldsol)$ % oldsol as vanishing expression
  164. prev_depend:=storedepend(vari)$
  165. h2:=sub(u_=oldsol_ex,pde);
  166. if h2 neq 0 then <<
  167. h1:=solve(oldsol_ex,vari);
  168. if h1 neq {} then <<
  169. h1=first h1;
  170. for each h3 in h1 do
  171. if freeof(h3,arbcomplex) then h2:=sub(h3,h2)
  172. >>
  173. >>;
  174. restoredepend(prev_depend)$
  175. if 0=h2 then return lisp
  176. <<write"The special solution to be generalized is an invariant ",
  177. "with respect to";terpri()$
  178. write"this symmetry, therefore no generalization is possible.";
  179. terpri()$
  180. for each h1 in fargs u_ do nodepend u_,h1;
  181. algebraic oldsol
  182. >>;
  183. pde:=pde-1;
  184. h1:= quasilinpde1(pde,u_,vari)$
  185. if h1 neq {} then <<
  186. h1:=first h1;
  187. % h2 is expressing the constants in terms of u_ and the xlist,ylist
  188. clist:={};
  189. h2:= for each e1 in h1 collect
  190. <<h3:=lisp(newfct(fname_,nil,nfct_))$
  191. lisp(nfct_:=add1 nfct_)$
  192. clist:=cons(h3,clist);
  193. h3 = e1>>;
  194. h20:=sub(u_=0,h2);
  195. h3:=solve(h2,vari);
  196. if h3 neq {} then <<
  197. h3:=first h3;
  198. h30:=sub(h20,h3);
  199. write"The substitutions to generalize the solution are: "$
  200. for each h4 in h30 do write h4$
  201. newsol:=sub(h30,oldsol_ex);
  202. % newsol:=second dropredundant(newsol,
  203. % alle Konstanten und
  204. % Funktionen die nicht in vari sind,
  205. % vari)$
  206. lisp <<write"The new solution";
  207. % if length algebraic h >2 then write"s are:" else write" is:";
  208. % terpri()$
  209. >>;
  210. write"0 = ",newsol
  211. >>
  212. >>;
  213. for each h1 in fargs u_ do nodepend u_,h1;
  214. return newsol
  215. end$ % of NewParam
  216. %----------------------------
  217. symbolic operator einfachst$
  218. symbolic procedure einfachst(a,x)$
  219. % a is an algebraic list and the element where x appears, but
  220. % appears simplest, is found
  221. begin
  222. scalar el1,el2,hp;
  223. hp:=10000;
  224. a:=cdr a;
  225. while a do <<
  226. el2:=car a;a:=cdr a;
  227. if not freeof(el2,x) and
  228. ((not el1 ) or
  229. (el2 = x ) or
  230. <<if not polyp(el2,cons(x,nil)) then nil
  231. else
  232. <<coeff1(el2,x,nil)$
  233. if hipow!*<hp then <<hp:=hipow!*; t>>
  234. else nil >> >> ) then el1:=el2
  235. >>;
  236. return el1
  237. end$
  238. %----------------------------
  239. algebraic procedure TransDf(y,yslist,vlist,indxlist)$
  240. begin
  241. scalar m,n,e1,dfy;
  242. return
  243. if indxlist={} then sub(yslist,y)
  244. else <<
  245. m:=first indxlist;
  246. n:=0;
  247. dfy:=TransDf(y,yslist,vlist,rest indxlist);
  248. for each e1 in vlist sum <<
  249. n:=n+1;
  250. df(dfy,e1)*Dv!/Dx(n,m)
  251. >>
  252. >>
  253. end$ % of TransDf
  254. %----------------------------
  255. algebraic procedure TransDeriv(yik,yslist,vlist)$
  256. begin
  257. scalar indxlist,y,l1,l2;
  258. indxlist:=lisp cons('LIST,combidif(yik))$
  259. return TransDf(first indxlist,yslist,vlist,rest indxlist)
  260. end$ % of TransDeriv
  261. %----------------------------
  262. algebraic procedure DeTrafo(eqlist,yslist,xslist,ulist,vlist)$
  263. % Transformations of all orders are performed (point-, contact-,...)
  264. % but only x-derivatives of y's are transformed. To transform other
  265. % any other derivatives, subdiff1 must be extended to include all other
  266. % occuring x-derivatives.
  267. begin
  268. scalar avar,nvar,detpd,n,m,ordr,e1,e2,e3,sb;
  269. m:=length(xslist); n:=length(yslist)+m;
  270. clear dyx!/duv,Dv!/Dx;
  271. matrix dyx!/duv(n,n);
  272. matrix Dv!/Dx(m,m);
  273. avar:=append(yslist,xslist);
  274. nvar:=append(ulist,vlist);
  275. n:=0;
  276. for each e1 in avar do <<
  277. n:=n+1;m:=0;
  278. for each e2 in nvar do <<
  279. m:=m+1;
  280. dyx!/duv(m,n):=df(rhs e1,e2)
  281. >>
  282. >>;
  283. detpd:=det(dyx!/duv);
  284. %write"detpd=",detpd;
  285. if detpd=0 then return
  286. <<write"The proposed transformation is not regular!";{}>>;
  287. clear dyx!/duv;
  288. ordr:=0;
  289. for each e1 in eqlist do
  290. for each e2 in yslist do
  291. <<n:=totdeg(e1,lhs e2);
  292. if n>ordr then ordr:=n>>;
  293. sb:=subdif1(for each e1 in xslist collect lhs e1,
  294. for each e1 in yslist collect lhs e1,ordr);
  295. % computation of Dv/Dx:=(Dx/Dv)^(-1)
  296. n:=0;
  297. for each e1 in xslist do <<
  298. n:=n+1;m:=0;
  299. for each e2 in vlist do <<
  300. m:=m+1;
  301. Dv!/Dx(n,m):=total_alg_mode_deriv(rhs e1,e2)
  302. % it is assumed ulist does depend on vlist
  303. >>
  304. >>;
  305. Dv!/Dx:=Dv!/Dx**(-1);
  306. % Substitution of all derivatives
  307. for each e1 in sb do
  308. for each e2 in e1 do <<
  309. if not freeof(eqlist,lhs e2) then <<
  310. % which function is to be differentiated wrt. which variable
  311. eqlist:=sub(lhs e2=TransDeriv(rhs e2,yslist,vlist),eqlist)
  312. >>
  313. >>;
  314. clear Dv!/Dx;
  315. return sub(xslist,sub(yslist,eqlist))$
  316. end$ % of DeTrafo
  317. %----------------------------
  318. algebraic procedure grouping(el1,el2,xlist,ylist,nx,ny)$
  319. begin scalar h,el3,xslist,yslist$
  320. %------- Grouping the new variables to ulist and vlist
  321. h:={};
  322. xslist:={}; % list of expressions to calculate new indep. var.
  323. yslist:={}; % list of expressions to calculate new dep. var.
  324. %---- at first the obvious allocations
  325. for each el3 in el1 do %-- all similarity variables
  326. if freeoflist(el3,ylist) then xslist:=cons(el3,xslist) else
  327. if freeoflist(el3,xlist) then yslist:=cons(el3,yslist) else
  328. h:=cons(el3,h);
  329. %---- now the symmetry variable
  330. if freeoflist(el2,ylist) or (length(yslist) = ny) then
  331. xslist:=cons(el2,xslist) else
  332. if freeoflist(el2,xlist) or (length(xslist) = nx) then
  333. yslist:=cons(el2,yslist) else
  334. xslist:=cons(el2,xslist);
  335. %---- now the remaining cases
  336. for each el3 in h do
  337. if length(yslist) < ny then yslist:=cons(el3,yslist)
  338. else xslist:=cons(el3,xslist);
  339. return {xslist,yslist}
  340. end$ % of grouping
  341. %----------------------------
  342. algebraic procedure rename_u_(xslist,yslist,el2,u_,u,v)$
  343. begin scalar i,vlist,ulist,el3,h,smv$
  344. %---- Renaming the u_ to ui in yslist and to vi in xslist
  345. i:=0;
  346. vlist:={};
  347. xslist:=for each el3 in xslist collect
  348. <<i:=i+1;
  349. if length xslist>1 then h:=mkid(v,i)
  350. else h:=v;
  351. vlist:=cons(h,vlist);
  352. if el3=el2 then smv:=h;
  353. sub(u_=h,el3)
  354. >>;
  355. i:=0;
  356. ulist:={};
  357. yslist:=for each el3 in yslist collect
  358. <<i:=i+1;
  359. if length yslist>1 then h:=mkid(u,i)
  360. else h:=u;
  361. ulist:=cons(h,ulist);
  362. if el3=el2 then smv:=h;
  363. sub(u_=h,el3)
  364. >>;
  365. return {xslist,yslist,reverse vlist,reverse ulist,smv}
  366. end$ % of rename_u_
  367. %----------------------------
  368. algebraic procedure solve_for_old_var(xslist,yslist,xlist,ylist,nx,ny)$
  369. begin scalar h1,h2$
  370. %---- Solve for old variables
  371. h1:=nil;
  372. h2:=solve(append(yslist,xslist),append(xlist,ylist));
  373. if h2={} then h1:=t
  374. else h2:=first h2; %--- possibly other solutions
  375. if LIST neq lisp(car algebraic h2) then h1:=t else
  376. if length(h2)<(nx+ny) then el2:=t;
  377. if h1 then repeat lisp
  378. <<write"The algebraic system ",append(xslist,yslist),
  379. " could not be solved for ",append(xlist,ylist),".";
  380. write"Please enter the solution in form of a list {",
  381. reval algebraic first xlist,"=...,...",
  382. reval algebraic first ylist,"=...,...} or enter a ",
  383. "semicolon ; to end this investigation:";
  384. algebraic(h2:=termxread())
  385. >> until h2=nil or ( lisp(pairp algebraic h2) and
  386. (LIST=lisp(car algebraic h2)) and
  387. (length(h2)=(nx+ny)) )
  388. else
  389. <<lisp<<terpri()$
  390. write"The suggested solution of the algebraic system which will";
  391. terpri()$
  392. write"do the transformation is: ";
  393. terpri()
  394. >>;
  395. write h2;
  396. if yesp "Is the solution ok?" then else lisp <<
  397. write"Please enter the solution in form of a list {",
  398. reval algebraic first xlist,"=...,...",
  399. reval algebraic first ylist,"=...,...} or enter a ",
  400. "semicolon ; to end this investigation:";
  401. algebraic(h2:=termxread())
  402. >>
  403. >>;
  404. return h2
  405. end$ % of solve_for_old_var$
  406. %----------------------------
  407. algebraic procedure switch_r_s(h2,smv,ylist,u,v)$
  408. begin scalar xslist,yslist,el3,h$
  409. %---- Exchange of dependent and independent variables
  410. xslist:={};
  411. yslist:={};
  412. for each el3 in h2 do if freeof(ylist,lhs el3) then
  413. xslist:=cons(el3,xslist) else
  414. yslist:=cons(el3,yslist);
  415. lisp <<terpri()$
  416. write"In the intended transformation shown above",
  417. " the dependent ";terpri()$
  418. if length yslist>2 then
  419. write"variables are the ",reval algebraic u,"i and " else
  420. write"variable is ",reval algebraic u," and ";
  421. if length xslist>2 then
  422. write"the independent variables are the ",
  423. reval algebraic v,"i." else
  424. write"the independent variable is ",reval algebraic v,".";
  425. terpri()$
  426. write"The symmetry variable is ",reval algebraic smv,", i.e. the ",
  427. "transformed expression";terpri();
  428. write"will be free of ",reval algebraic smv,".";
  429. >>;
  430. h:=if yesp "Is this selection of dependent and independent variables ok?"
  431. then nil else <<
  432. lisp <<write"Please enter a list of substitutions. For example, to";
  433. terpri()$
  434. write"make the variable, which is so far call u1, to an";
  435. terpri()$
  436. write"independent variable v2 and the variable, which is ";
  437. terpri()$
  438. write"so far called v2, to an dependent variable u1, ";
  439. terpri()$
  440. write"enter: `{u1=v2, v2=u1};'">>;
  441. termxread()
  442. >>;
  443. if h and (h neq {}) then <<xslist:=sub(h,xslist);
  444. yslist:=sub(h,yslist);
  445. smv :=sub(h,smv)>>;
  446. return {xslist,yslist,smv}
  447. end$ % of switch_r_s
  448. %----------------------------
  449. algebraic procedure makepde(genlist,u_)$
  450. begin scalar h,el2,el3,vari,bv;
  451. vari:={};
  452. return
  453. cons(
  454. num for each el2 in genlist sum
  455. <<h:=lhs el2;
  456. h:=lisp <<
  457. el3:=explode reval algebraic h;
  458. bv:=t;
  459. while bv do <<
  460. if car el3 ='!_ then bv:=nil;
  461. el3:=cdr el3
  462. >>;
  463. intern compress el3
  464. >>;
  465. depend u_,h;
  466. vari:=cons(h,vari);
  467. (rhs el2) * df(u_,h)
  468. >>,
  469. vari)
  470. end$ % of makepde
  471. %----------------------------
  472. algebraic procedure totdeglist(eqlist,ylist)$
  473. begin scalar n,ordr,e1,e2;
  474. ordr:=0;
  475. for each e1 in eqlist do
  476. for each e2 in ylist do
  477. <<n:=totdeg(e1,e2);
  478. if n>ordr then ordr:=n>>;
  479. return ordr
  480. end$ % of totdeglist
  481. %----------------------------
  482. algebraic procedure similarity(problem,genlist,con,u,v)$
  483. % con ... the free constants/functions in the general symmetry
  484. % u ... the name of the new independent variables
  485. % v ... the name of the new dependent variables
  486. begin scalar vari,pde,el1,el2,el3,el4,copgen,symvarfound,
  487. trans1,trans2,i,j,h,h2,n,denew,xlist,ylist,
  488. eqlist,ulist,vlist,nx,ny,xslist,yslist,smv,
  489. trafoprob;
  490. cpu:=lisp time()$ gc:=lisp gctime()$
  491. %--------- extracting input data
  492. eqlist:=maklist first problem;
  493. ylist :=maklist second problem; ny:=length ylist;
  494. xlist :=maklist third problem; nx:=length xlist;
  495. trafoprob:={problem,nil}; % to be returned if trafo not possible
  496. problem:=nil;
  497. eqlist:=for each el1 in eqlist collect equ_to_expr el1;
  498. % if length eqlist > 1 then eqlist:=desort eqlist;
  499. %--------- initial printout
  500. lisp(
  501. if tr_as then terpri());
  502. %--------- initializations
  503. ordr:=totdeglist(eqlist,ylist)$
  504. vari:=append(ylist,xlist);
  505. for each el1 in xlist do
  506. for each el2 in ylist do
  507. if not my_freeof(el2,el1) then nodepend el2,el1;
  508. lisp(
  509. if tr_as then <<
  510. write "The ODE/PDE (-system) under investigation is :";terpri()$
  511. for each el1 in cdr eqlist do algebraic write"0 = ",el1;
  512. terpri()$write "for the function(s) : ";
  513. fctprint( cdr reval algebraic ylist);write".";
  514. terpri()$terpri()
  515. >>);
  516. lisp(
  517. if tr_as then <<
  518. if length ylist >2 then % not >1 because of alg. list in symb. mode
  519. write"It will be looked for new dependent variables ",u,"i "
  520. else
  521. write"It will be looked for a new dependent variable ",u;
  522. terpri()$
  523. if length xlist >2 then
  524. write"and independent variables ",v,"i"
  525. else
  526. write"and an independent variable ",v;
  527. write" such that the transformed";
  528. terpri()$
  529. write"de(-system) does not depend on ",u;
  530. if length ylist >2 then write"1";
  531. write" or ",v$
  532. if length xlist >2 then write"1";
  533. write".";
  534. terpri()
  535. >>);
  536. % for each el1 in con do % i.e. for all symmetries do:
  537. % if freeoflist(el1,xlist) and
  538. % freeoflist(el1,ylist) then %------ no pseudo-Lie-symm.
  539. % <<copgen:=sub(el1=1,genlist); %------ calculate symmetry
  540. % for each el2 in con do
  541. % if el1 neq el2 then copgen:=sub(el2=0,copgen);
  542. copgen:=genlist;
  543. % write"The symmetry now under investigation is:";
  544. % for each el1 in copgen do write el1;
  545. %---------- formulate the PDE for the similarity variables
  546. pde:=first makepde(copgen,u_);
  547. %--------- find similarity variable
  548. trans2 :={};
  549. lisp<<terpri()$write"1. Determination of the similarity variable";
  550. if nx+ny>2 then write"s">>;
  551. trans1 := quasilinpde1(pde,u_,vari); % for the similarity variable
  552. if trans1 neq {} then
  553. <<
  554. %-------------- Determining the similarity variables ui_
  555. i:=0;
  556. trans1:=for each el1 in trans1 collect
  557. <<i:=i+1;
  558. h:=length(genlist)-1;
  559. if h=1 then % one single ODE
  560. <<
  561. % write"In the following 1 similarity variable U_ has to be";
  562. % write"determined through 0 = ff where ff is an"$
  563. % write"arbitrary function of arguments given in the ",
  564. % "following list:"$
  565. % write el1;
  566. el2:=num(first el1 - second el1);
  567. if freeof(el2,u_) then
  568. el2:=num(first el1 - 2*second el1);
  569. lisp<<write"A suggestion for this function ff provides:"$
  570. terpri()>>$
  571. write"0 = ", el2$
  572. if yesp "Do you like this choice?" then {el2}
  573. else <<
  574. repeat <<
  575. lisp <<
  576. write"Put in an alternative expression which "$terpri()$
  577. write"- is functionally dependent only on elements of",
  578. " ff given above and"$ terpri()$
  579. write"- depends on U_ and if set to zero determines U_"$
  580. terpri()>>$
  581. h:=termxread()$
  582. >> until not freeof(h,U_)$
  583. {h}
  584. >>
  585. >> else % a PDE or a system of DEs
  586. lisp
  587. <<
  588. % terpri()$
  589. % write"Now the similarity variables U_i, i=1,...",h,
  590. % " have to be"$terpri()$
  591. % write"determined through conditions 0 = ffi, i=1,...",h,
  592. % ", where ffi are"$terpri()$
  593. % write"arbitrary functions ";terpri()$
  594. % algebraic(write "ffi = ",
  595. % lisp( cons('ffi,cdr reval algebraic el1)));
  596. % terpri()$
  597. % write"such that the functional determinant of these ",
  598. % "expressions"$ terpri()$
  599. % write"including u_ from above taken w.r.t. ";
  600. % for each el3 in cdr algebraic ylist do write reval el3,",";
  601. % write reval cadr algebraic xlist;
  602. % for each el3 in cddr algebraic xlist do write ",",reval el3;
  603. % terpri()$
  604. % write"must not vanish."$terpri()$
  605. algebraic <<
  606. el2:=einfachst(el1,u_);
  607. h2:={};
  608. for each el3 in el1 do
  609. if el3 neq el2 then h2:=cons(num(el2-el3),h2)
  610. >>;
  611. write"A suggestion for these functions ffi in form of a list ",
  612. "{ff1,ff2,... } is: "$terpri()$
  613. deprint(cdr reval algebraic h2);
  614. if yesp "Do you like this choice?" then algebraic h2
  615. else <<
  616. write"Put in an alternative list of expression which"$
  617. terpri()$
  618. write"- are functionally dependent only on the above ",
  619. "arguments and"$terpri()$
  620. write"- which if set to zero determine U_i, i.e."$terpri()$
  621. write"- the functional determinant of these expressions"$
  622. terpri()$
  623. write" including U_ from above taken w.r.t. ",
  624. cdr reval algebraic append(ylist,xlist)$
  625. terpri()$
  626. write" must not vanish."$terpri()$
  627. algebraic(h2):=termxread()
  628. >>
  629. >>
  630. >>$
  631. %--------- find symmetry variable
  632. pde:=pde-1;
  633. lisp<<terpri()$write"2. Determination of the symmetry variable">>;
  634. trans2 := quasilinpde1(pde,u_,vari); % for the symmetry variable
  635. if (length xlist=1) and (trans2={}) then
  636. for each e1 in fargs u_ do nodepend u_,e1
  637. else
  638. <<% If no symmetry variable is found (trans2={}) then proceed
  639. % only if special solutions of PDEs are to be found with the
  640. % solution being only a function of the similarity variables
  641. if trans2={} then << % take any variable
  642. h:=reverse vari; % reverse to have not a function as symvar
  643. while 1+sub(u_=first h,pde)=0 do h:=rest h; % no similarity var.
  644. for each e1 in fargs u_ do nodepend u_,e1;
  645. h:=first h;
  646. trans2:={u_ - h};
  647. lisp<<write"Because the correct symmetry variable was not ",
  648. "found, the program will";terpri()$
  649. write"take ",reval algebraic h,
  650. " instead with the consequence ",
  651. "that not the whole transformed ";terpri()$
  652. write"PDE will be free of ",
  653. reval algebraic h," but only those ",
  654. "terms without ",
  655. reval algebraic h,"-derivative";terpri()$
  656. write"which is still of use for finding special ",
  657. reval algebraic h,"-independent solutions ";
  658. terpri()$
  659. write"of the PDE."
  660. >>
  661. >> else <<
  662. %--------------- Determining an optimal symmetry variable
  663. for each e1 in fargs u_ do nodepend u_,e1;
  664. symvarfound:=t;
  665. i:=0;
  666. trans2:=for each el1 in trans2 collect
  667. <<i:=i+1;
  668. lisp
  669. <<
  670. % terpri()$
  671. % write"In the following the symmetry variable U_",
  672. % " has to be"$terpri()$
  673. % write"determined through a condition 0 = ff",
  674. % ", where ff is"$terpri()$
  675. % write"an arbitrary function ";terpri()$
  676. % algebraic(write"ff = ",
  677. % lisp( cons('ff,cdr reval algebraic el1)));terpri()$
  678. write"A suggestion for this function ff(..) yields:";terpri()
  679. >>;
  680. h:=einfachst(el1,u_);
  681. if lisp<<h2:=reval algebraic(num h);
  682. (not pairp h2) or (car h2 neq 'PLUS)>> then
  683. h:=num(h+1);
  684. % if h= first el1 then
  685. % if freeof(num(h+second el1),u_) then h:=num(h+2*second el1)
  686. % else h:=num(h+ second el1)
  687. % else
  688. % if freeof(num(h+ first el1),u_) then h:=num(h+2* first el1)
  689. % else h:=num(h+ first el1);
  690. write"0 = ",h$
  691. if yesp "Do you like this choice?" then h
  692. else <<
  693. repeat <<
  694. lisp <<
  695. write"Put in an alternative expression which "$terpri()$
  696. write"- is functionally dependent only on arguments of",
  697. " ff given above and"$terpri()$
  698. write"- depends on u_ and if set to zero determines u_"$
  699. terpri()>>$
  700. h:=termxread()$
  701. >> until not freeof(h,u_)$
  702. h
  703. >>
  704. >>
  705. >>$
  706. %for each el1 in trans1 do
  707. %for each el2 in trans2 do
  708. el1:=first trans1;
  709. el2:=first trans2;
  710. % <<
  711. %------- Grouping the new variables to ulist and vlist
  712. yslist:=grouping(el1,el2,xlist,ylist,nx,ny)$
  713. xslist:= first yslist;
  714. yslist:=second yslist;
  715. %---- Renaming the u_ to ui in yslist and to vi in xslist
  716. smv:=rename_u_(xslist,yslist,el2,u_,u,v)$
  717. xslist:= first smv;
  718. yslist:=second smv;
  719. vlist := third smv;
  720. ulist:=part(smv,4);
  721. smv:=part(smv,5); % the symmetry variable
  722. %---- Solve for old variables
  723. h2:=solve_for_old_var(xslist,yslist,xlist,ylist,nx,ny);
  724. if h2 neq nil then <<
  725. %---- Exchange of dependent and independent variables
  726. smv:=switch_r_s(h2,smv,ylist,u,v)$
  727. xslist:= first smv;
  728. yslist:=second smv;
  729. smv :=third smv;
  730. %---- Doing the point transformation
  731. for each el3 in ulist do
  732. <<for each el4 in fargs el3 do nodepend el3,el4;
  733. for each el4 in vlist do
  734. %if el4 neq smv then % if new DEs without symm. var. smv
  735. depend el3,el4>>;
  736. for each el3 in ylist do
  737. for each el4 in xlist do depend el3,el4;
  738. eqlist:=DeTrafo(eqlist,yslist,xslist,ulist,vlist);
  739. lisp(
  740. if tr_as then <<
  741. terpri()$
  742. write"The transformed equation";
  743. if length(algebraic eqlist)>2 then write"s";
  744. if symvarfound then
  745. write" which should be free of ",reval algebraic smv,":"
  746. else
  747. write" in which the terms without ",reval algebraic smv,
  748. "-derivative are free of ",reval algebraic smv,":";
  749. terpri()
  750. >>);
  751. eqlist:=for each el3 in eqlist collect <<
  752. el3:=factorize num el3;
  753. for each el4 in el3 product
  754. if 0=totdeglist({el4},ulist) then 1
  755. else el4
  756. >>;
  757. lisp deprint(cdr reval algebraic eqlist);
  758. if (length(vlist)>1) and (not freeof(vlist,smv)) then <<
  759. vlist:=cons(smv,lisp(delete(reval algebraic smv,
  760. reval algebraic vlist)));
  761. if yesp
  762. "Shall the dependence on the symmetry variable be dropped?"
  763. then
  764. <<for each el3 in ulist do
  765. if not my_freeof(el3,smv) then nodepend el3,smv;
  766. vlist:=rest vlist>>;
  767. eqlist:=for each el3 in eqlist collect <<
  768. el3:=factorize num el3;
  769. for each el4 in el3 product
  770. if 0=totdeglist({el4},ulist) then 1
  771. else el4
  772. >>
  773. >>;
  774. trafoprob:={{eqlist,ulist,vlist},append(xslist,yslist)}
  775. >>
  776. % >>
  777. >>
  778. >>;
  779. % >>;
  780. for each el1 in xlist do
  781. for each el2 in ylist do
  782. depend el2,el1;
  783. clear ff,ffi;
  784. return trafoprob
  785. end$ % of similarity
  786. %----------------------------
  787. algebraic procedure quasilinpde1(pde,u_,vari)$
  788. begin scalar trans1,e1,e2,q;
  789. trans1 := quasilinpde(pde,u_,vari); % for the similarity variable
  790. if trans1={} then <<
  791. write"The program was not able to find the general solution ",
  792. "of the PDE: ",pde," for the function ",u_,".";
  793. lisp <<
  794. write"Please enter either only a semicolon if no solution ",
  795. "is known or enter "$terpri()$
  796. write"the solution of the PDE in form ",
  797. "of a list {A1,A2,...} where ";terpri()$
  798. write"the Ai are algebraic expressions in ",
  799. cdr reval algebraic cons(u_,vari);terpri()$
  800. write"such that any function ff(A1,A2,...) which is not ",
  801. "independent of `",u_;write"'";terpri()$
  802. write"determines a solution `",u_,"' of the PDE through 0=ff: "
  803. >>;
  804. trans1:={termxread()};
  805. if trans1={nil} then trans1:={}
  806. >>;
  807. return trans1
  808. end$ % of quasilinpde1
  809. end$