# source:git/Singular/LIB/findifs.lib@3686937

spielwiese
Last change on this file since 3686937 was 3686937, checked in by Oleksandr Motsak <motsak@…>, 10 years ago
Added '$Id$' as a comment to all libs (LIB/*.lib)
• Property mode set to 100644
File size: 18.1 KB
Line
1//////////////////////////////////////////////////////////////////////////////
2version="version findifs.lib 4.0.0.0 Jun_2013 "; // $Id$
3category="Applications";
4info="
5LIBRARY: findifs.lib  Tools for the finite difference schemes
6AUTHORS: Viktor Levandovskyy,     levandov@math.rwth-aachen.de
7
8OVERVIEW:
9        We provide the presentation of difference operators in a polynomial,
10        semi-factorized and a nodal form. Running @code{findifs_example();}
11        will demonstrate, how we generate finite difference schemes of linear PDEs
12        from given approximations.
13
14Theory: The method we use have been developed by V. Levandovskyy and Bernd Martin. The
15computation of a finite difference scheme of a given single linear partial
16differential equation with constant coefficients with a given approximation
17rules boils down to the computation of a Groebner basis of a submodule of
18a free module with respect to the ordering, eliminating module components.
19
20
21Support: SpezialForschungsBereich F1301 of the Austrian FWF
22
23PROCEDURES:
24findifs_example();  containes a guided explanation of our approach
25decoef(P,n);  decompose polynomial P into summands with respect to the number n
26difpoly2tex(S,P[,Q]); present the difference scheme in the nodal form
27
28
29exp2pt(P[,L]); convert a polynomial M into the TeX format, in nodal form
30texcoef(n);  converts the number n into TeX
31npar(n); search for 'n' among the parameters and returns its number
32magnitude(P); compute the square of the magnitude of a complex expression
33replace(s,what,with);  replace in s all the substrings with a given string
34xchange(w,a,b); exchange two substrings in a given string
35
37";
38
39
40LIB "latex.lib";
41LIB "poly.lib";
42
43proc tst_findif()
44{
45  example decoef;
46  example difpoly2tex;
47  example exp2pt;
48  example texcoef;
49  example npar;
50  example magnitude;
51  example replace;
52  example xchange;
53}
54
55// static procs:
56// par2tex(s);     convert special characters to TeX in s
57// mon2pt(P[,L]); convert a monomial M into the TeX format, in nodal form
58
59
60// 1. GLOBAL ASSUME: in the ring we have first Tx, then Tt: [FIXED, not needed anymore]!
61// 2. map vars other than Tx,Tt to parameters instead or just ignore them [?]
62// 3. clear the things with brackets
63// 4. todo: content resp lcmZ, gcdZ
64
65proc xchange(string where, string a, string b)
66"USAGE:  xchange(w,a,b);  w,a,b strings
67RETURN:  string
68PURPOSE:  exchanges substring 'a' with a substring 'b' in the string w
69NOTE:
70EXAMPLE: example xchange; shows examples
71"{
72  // replaces a<->b in where
73  // assume they are of the same size [? seems to work]
74  string s = "H";
75  string t;
76  t = replace(where,a,s);
77  t = replace(t,b,a);
78  t = replace(t,s,b);
79  return(t);
80}
81example
82{
83 " EXAMPLE:"; echo=2;
84 ring r = (0,dt,dh,A),Tt,dp;
85 poly p = (Tt*dt+dh+1)^2+2*A;
86 string s = texpoly("",p);
87 s;
88 string t = xchange(s,"dh","dt");
89 t;
90}
91
92static proc par2tex(string s)
93"USAGE:  par2tex(s);  s a string
94RETURN:  string
95PURPOSE: converts special characters to TeX in s
96NOTE: the convention is the following:
97      'Tx' goes to 'T_x',
98      'dx' to '\\tri x' (the same for dt, dy, dz),
99      'theta', 'ro', 'A', 'V' are converted to greek letters.
100EXAMPLE: example par2tex; shows examples
101"{
102  // can be done with the help of latex_lib
103
104// s is a tex string with a poly
105// replace theta with \theta
106// A with \lambda
107// dt with \tri t
108// dh with \tri h
109// Tx with T_x, Ty with T_y
110// Tt with T_t
111// V with \nu
112// ro with \rho
113// dx with \tri x
114// dy with \tri y
115  string t = s;
116  t = replace(t,"Tt","T_t");
117  t = replace(t,"Tx","T_x");
118  t = replace(t,"Ty","T_y");
119  t = replace(t,"dt","\\tri t");
120  t = replace(t,"dh","\\tri h");
121  t = replace(t,"dx","\\tri x");
122  t = replace(t,"dy","\\tri y");
123  t = replace(t,"theta","\\theta");
124  t = replace(t,"A","\\lambda");
125  t = replace(t,"V","\\nu");
126  t = replace(t,"ro","\\rho");
127  return(t);
128}
129example
130{
131 " EXAMPLE:"; echo=2;
132 ring r = (0,dt,theta,A),Tt,dp;
133 poly p = (Tt*dt+theta+1)^2+2*A;
134 string s = texfactorize("",p);
135 s;
136 par2tex(s);
137 string T = texfactorize("",p*(-theta*A));
138 par2tex(T);
139}
140
141proc replace(string s, string what, string with)
142"USAGE:  replace(s,what,with);  s,what,with strings
143RETURN:  string
144PURPOSE: replaces in 's' all the substrings 'what' with substring 'with'
145NOTE:
146EXAMPLE: example replace; shows examples
147"{
148  // clear: replace in s, "what" with "with"
149  int ss = size(s);
150  int cn = find(s,what);
151  if ( (cn==0) || (cn>ss))
152  {
153    return(s);
154  }
155  int gn = 0; // global counter
156  int sw = size(what);
157  int swith = size(with);
158  string out="";
159  string tmp;
160  gn  = 0;
161  while(cn!=0)
162  {
163//    "cn:";    cn;
164//    "gn";    gn;
165    tmp = "";
166    if (cn>gn)
167    {
168      tmp = s[gn..cn-1];
169    }
170//    "tmp:";tmp;
171//    out = out+tmp+" "+with;
172    out = out+tmp+with;
173//    "out:";out;
174    gn  = cn + sw;
175    if (gn>ss)
176    {
177      // ( (gn>ss) || ((sw>1) && (gn >= ss)) )
178      // no need to append smth
179      return(out);
180    }
181//     if (gn == ss)
182//     {
183
184//     }
185    cn  = find(s,what,gn);
186  }
187  // and now, append the rest of s
188  //  out = out + " "+ s[gn..ss];
189  out = out + s[gn..ss];
190  return(out);
191}
192example
193{
194 " EXAMPLE:"; echo=2;
195 ring r = (0,dt,theta),Tt,dp;
196 poly p = (Tt*dt+theta+1)^2+2;
197 string s = texfactorize("",p);
198 s;
199 s = replace(s,"Tt","T_t"); s;
200 s = replace(s,"dt","\\tri t"); s;
201 s = replace(s,"theta","\\theta"); s;
202}
203
204proc exp2pt(poly P, list #)
205"USAGE:  exp2pt(P[,L]);  P poly, L an optional list of strings
206RETURN:  string
207PURPOSE: convert a polynomial M into the TeX format, in nodal form
208ASSUME: coefficients must not be fractional
209NOTE: an optional list L contains a string, which will replace the default
210value 'u' for the discretized function
211EXAMPLE: example exp2pt; shows examples
212"{
213// given poly in vars [now Tx,Tt are fixed],
214// create Tex expression for points of lattice
215// coeffs must not be fractional
216  string varnm = "u";
217  if (size(#) > 0)
218  {
219    if (typeof(#[1])=="string")
220    {
221      varnm = string(#[1]);
222    }
223  }
224  //  varnm;
225  string rz,mz;
226  while (P!=0)
227  {
228    mz = mon2pt(P,varnm);
229    if (mz[1]=="-")
230    {
231      rz = rz+mz;
232    }
233    else
234    {
235      rz = rz + "+" + mz;
236    }
238  }
239  rz  = rz[2..size(rz)];
240  return(rz);
241}
242example
243{
244  " EXAMPLE:"; echo=2;
245  ring r = (0,dh,dt),(Tx,Tt),dp;
246  poly M = (4*dh*Tx^2+1)*(Tt-1)^2;
247  print(exp2pt(M));
248  print(exp2pt(M,"F"));
249}
250
251static proc mon2pt(poly M, string V)
252"USAGE:  mon2pt(M,V);  M poly, V a string
253RETURN:  string
254PURPOSE: convert a monomial M into the TeX format, nodal form
255EXAMPLE: example mon2pt; shows examples
256"{
257  // searches for Tx, then Tt
258  // monomial to the lattice point conversion
259  // c*X^a*Y^b --> c*U^{n+a}_{j+b}
262  //  int a = e[2]; // convention: first Tx, then Tt
263  //  int b = e[1];
264  int i;
265  int a , b, c = 0,0,0;
266  int ia,ib,ic = 0,0,0;
267  int nv = nvars(basering);
268  string s;
269  for (i=1; i<=nv ; i++)
270  {
271    s = string(var(i));
272    if (s=="Tt") { a = e[i]; ia = i;}
273    if (s=="Tx") { b = e[i]; ib = i;}
274    if (s=="Ty") { c = e[i]; ic = i;}
275  }
279  //  string tc = texobj("",c); // why not texpoly?
280  string tc = texcoef(cM);
281  string rs;
282  if (cM==-1)
283  {
284    rs = "-";
285  }
286  if (cM^2 != 1)
287  {
288// we don't need 1 or -1 as coeffs
289//    rs = clTex(tc)+" ";
290//    rs = par2tex(rmDol(tc))+" ";
291    rs = par2tex(tc)+" ";
292  }
293// a = 0 or b = 0
294  rs = rs + V +"^{n";
295  if (a!=0)
296  {
297    rs = rs +"+"+string(a);
298  }
299  rs = rs +"}_{j";
300  if (b!=0)
301  {
302    rs = rs +"+"+string(b);
303  }
304  if (c!=0)
305  {
306    rs = rs + ",k+";
307    rs = rs + string(c);
308  }
309  rs = rs +"}";
310  return(rs);
311}
312example
313{
314 "EXAMPLE:"; echo=2;
315  ring r = (0,dh,dt),(Tx,Tt),dp;
316  poly M = (4*dh^2-dt)*Tx^3*Tt;
317  print(mon2pt(M,"u"));
318  poly N = ((dh-dt)/(dh+dt))*Tx^2*Tt^2;
319  print(mon2pt(N,"f"));
320  ring r2 = (0,dh,dt),(Tx,Ty,Tt),dp;
321  poly M  = (4*dh^2-dt)*Tx^3*Ty^2*Tt;
322  print(mon2pt(M,"u"));
323}
324
325proc npar(number n)
326"USAGE:  npar(n);  n a number
327RETURN:  int
328PURPOSE: searches for 'n' among the parameters and returns its number
329EXAMPLE: example npar; shows examples
330"{
331  // searches for n amongst parameters
332  // and returns its number
333  int i,j=0,0;
334  list L = ringlist(basering);
335  list M = L[1][2]; // pars
336  string sn = string(n);
337  sn = sn[2..size(sn)-1];
338  for (i=1; i<=size(M);i++)
339  {
340    if (M[i] == sn)
341    {
342      j = i;
343    }
344  }
345  if (j==0)
346  {
347    "Incorrect parameter";
348  }
349  return(j);
350}
351example
352{
353 "EXAMPLE:"; echo=2;
354  ring r = (0,dh,dt,theta,A),t,dp;
355  npar(dh);
356  number T = theta;
357  npar(T);
358  npar(dh^2);
359}
360
361proc decoef(poly P, number n)
362"USAGE:  decoef(P,n);  P a poly, n a number
363RETURN:  ideal
364PURPOSE:  decompose polynomial P into summands with respect
365to the  presence of the number n in the coefficients
366NOTE:    n is usually a parameter with no power
367EXAMPLE: example decoef; shows examples
368"{
369  // decomposes polynomial into summands
370  // w.r.t. the presence of a number n in coeffs
371  // returns ideal
372  def br = basering;
373  int i,j=0,0;
374  int pos = npar(n);
375  if ((pos==0) || (P==0))
376  {
377    return(0);
378  }
379  pos = pos + nvars(basering);
380// map all pars except to vars, provided no things are in denominator
381  number con = content(P);
382  con = numerator(con);
383  P = cleardenom(P); //destroys content!
384  P = con*P; // restore the numerator part of the content
385  list M = ringlist(basering);
386  list L = M[1..4];
387  list Pars = L[1][2];
388  list Vars = L[2] + Pars;
389  L[1] = L[1][1]; // characteristic
390  L[2] = Vars;
391  // for non-comm things: don't need nc but graded algebra
392  //  list templ;
393  //  L[5] = templ;
394  //  L[6] = templ;
395  def @R = ring(L);
396  setring @R;
397  poly P = imap(br,P);
398  poly P0 = subst(P,var(pos),0);
399  poly P1 = P - P0;
400  ideal I = P0,P1;
401  setring br;
402  ideal I = imap(@R,I);
403  kill @R;
404  // check: P0+P1==P
405  poly Q = I[1]+I[2];
406  if (P!=Q)
407  {
408    "Warning: problem in decoef";
409  }
410  return(I);
411// substract the pure part from orig and check if n is remained there
412}
413example
414{
415 " EXAMPLE:"; echo=2;
416  ring r = (0,dh,dt),(Tx,Tt),dp;
417  poly P = (4*dh^2-dt)*Tx^3*Tt + dt*dh*Tt^2 + dh*Tt;
418  decoef(P,dt);
419  decoef(P,dh);
420}
421
422proc texcoef(number n)
423"USAGE:  texcoef(n);  n a number
424RETURN:  string
425PURPOSE:  converts the number n into TeX format
426NOTE: if n is a polynomial, texcoef adds extra brackets and performs some space substitutions
427EXAMPLE: example texcoef; shows examples
428"{
429  // makes tex from n
430  // and uses substitutions
431  // if n is a polynomial, adds brackets
432  number D  = denominator(n);
433  int DenIsOne = 0;
434  if ( D==number(1) )
435  {
436    DenIsOne = 1;
437  }
438  string sd = texpoly("",D);
439  sd = rmDol(sd);
440  sd = par2tex(sd);
441  number N  = numerator(n);
442  string sn = texpoly("",N);
443  sn = rmDol(sn);
444  sn = par2tex(sn);
445  string sout="";
446  int i;
447  int NisPoly = 0;
448  if (DenIsOne)
449  {
450    sout = sn;
451    for(i=1; i<=size(sout); i++)
452    {
453      if ( (sout[i]=="+") || (sout[i]=="-") )
454      {
455        NisPoly = 1;
456      }
457    }
458    if (NisPoly)
459    {
460      sout = "("+sout+")";
461    }
462  }
463  else
464  {
465    sout = "\\frac{"+sn+"}{"+sd+"}";
466  }
467  return(sout);
468}
469example
470{
471 " EXAMPLE:"; echo=2;
472  ring r = (0,dh,dt),(Tx,Tt),dp;
473  number n1,n2,n3 = dt/(4*dh^2-dt),(dt+dh)^2, 1/dh;
474  n1; texcoef(n1);
475  n2; texcoef(n2);
476  n3; texcoef(n3);
477}
478
479static proc rmDol(string s)
480{
481  // removes $s and _no_ (s on appearance 482 int i = size(s); 483 if (s[1] == "$") { s = s[2..i]; i--;}
484  if (s[1] == "(") { s = s[2..i]; i--;}
485  if (s[i] == "$") { s = s[1..i-1]; i--;} 486 if (s[i] == ")") { s = s[1..i-1];} 487 return(s); 488} 489 490proc difpoly2tex(ideal S, list P, list #) 491"USAGE: difpoly2tex(S,P[,Q]); S an ideal, P and optional Q are lists 492RETURN: string 493PURPOSE: present the difference scheme in the nodal form 494ASSUME: ideal S is the result of @code{decoef} procedure 495NOTE: a list P may be empty or may contain parameters, which will not 496appear in denominators 497@* an optional list Q represents the part of the scheme, depending 498on other function, than the major part 499EXAMPLE: example difpoly2tex; shows examples 500" 501{ 502 // S = sum s_i = orig diff poly or 503 // the result of decoef 504 // P = list of pars (numbers) not to be divided with, may be empty 505 // # is an optional list of polys, repr. the part dep. on "f", not on "u" 506 // S = simplify(S,2); // destroys the leadcoef 507 // rescan S and remove 0s from it 508 int i; 509 ideal T; 510 int ss = ncols(S); 511 int j=1; 512 for(i=1; i<=ss; i++) 513 { 514 if (S[i]!=0) 515 { 516 T[j]=S[i]; 517 j++; 518 } 519 } 520 S = T; 521 ss = j-1; 522 int GotF = 1; 523 list F; 524 if (size(#)>0) 525 { 526 F = #; 527 if ( (size(F)==1) && (F[1]==0) ) 528 { 529 GotF = 0; 530 } 531 } 532 else 533 { 534 GotF = 0; 535 } 536 int sf = size(F); 537 538 ideal SC; 539 int sp = size(P); 540 intvec np; 541 int GotP = 1; 542 if (sp==0) 543 { 544 GotP = 0; 545 } 546 if (sp==1) 547 { 548 if (P[1]==0) 549 { 550 GotP = 0; 551 } 552 } 553 if (GotP) 554 { 555 for (i=1; i<=sp; i++) 556 { 557 np[i] = npar(P[i])+ nvars(basering); 558 } 559 } 560 for (i=1; i<=ss; i++) 561 { 562 SC[i] = leadcoef(S[i]); 563 } 564 if (GotF) 565 { 566 for (i=1; i<=sf; i++) 567 { 568 SC[ss+i] = leadcoef(F[i]); 569 } 570 } 571 def br = basering; 572// map all pars except to vars, provided no things are in denominator 573 list M = ringlist(basering); 574 list L = M[1..4]; // erase nc part 575 list Pars = L[1][2]; 576 list Vars = L[2] + Pars; 577 L[1] = L[1][1]; // characteristic 578 L[2] = Vars; 579 580 def @R = ring(L); 581 setring @R; 582 ideal SC = imap(br,SC); 583 if (GotP) 584 { 585 for (i=1; i<=sp; i++) 586 { 587 SC = subst(SC,var(np[i]),1); 588 } 589 } 590 poly q=1; 591 q = lcm(q,SC); 592 setring br; 593 poly q = imap(@R,q); 594 number lq = leadcoef(q); 595 // lq; 596 number tmp; 597 string sout=""; 598 string vname = "u"; 599 for (i=1; i<=ss; i++) 600 { 601 tmp = leadcoef(S[i]); 602 S[i] = S[i]/tmp; 603 tmp = tmp/lq; 604 sout = sout +"+ "+texcoef(tmp)+"\\cdot ("+exp2pt(S[i])+")"; 605 } 606 if (GotF) 607 { 608 vname = "p"; //"f"; 609 for (i=1; i<=sf; i++) 610 { 611 tmp = leadcoef(F[i]); 612 F[i] = F[i]/tmp; 613 tmp = tmp/lq; 614 sout = sout +"+ "+texcoef(tmp)+"\\cdot ("+exp2pt(F[i],vname)+")"; 615 } 616 } 617 sout = sout[3..size(sout)]; //rm first + 618 return(sout); 619} 620example 621{ 622 "EXAMPLE:"; echo=2; 623 ring r = (0,dh,dt,V),(Tx,Tt),dp; 624 poly M = (4*dh*Tx+dt)^2*(Tt-1) + V*Tt*Tx; 625 ideal I = decoef(M,dt); 626 list L; L[1] = V; 627 difpoly2tex(I,L); 628 poly G = V*dh^2*(Tt-Tx)^2; 629 difpoly2tex(I,L,G); 630} 631 632 633proc magnitude(poly P) 634"USAGE: magnitude(P); P a poly 635RETURN: poly 636PURPOSE: compute the square of the magnitude of a complex expression 637ASSUME: i is the variable of a basering 638EXAMPLE: example magnitude; shows examples 639" 640{ 641 // check whether i is present among the vars 642 list L = ringlist(basering)[2]; // vars 643 int j; int cnt = 0; 644 for(j=size(L);j>0;j--) 645 { 646 if (L[j] == "i") 647 { 648 cnt = 1; break; 649 } 650 } 651 if (!cnt) 652 { 653 ERROR("a variable called i is expected in basering"); 654 } 655 // i is present, check that i^2+1=0; 656// if (NF(i^2+1,std(0)) != 0) 657// { 658// "Warning: i^2+1=0 does not hold. Reduce the output manually"; 659// } 660 poly re = subst(P,i,0); 661 poly im = (P - re)/i; 662 return(re^2+im^2); 663} 664example 665{ 666 "EXAMPLE:"; echo=2; 667 ring r = (0,d),(g,i,sin,cos),dp; 668 poly P = d*i*sin - g*cos +d^2*i; 669 NF( magnitude(P), std(i^2+1) ); 670} 671 672 673static proc clTex(string s) 674// removes beginning and ending$'s
675{
676  string t;
677  if (size(s)>2)
678  {
679    // why -3?
680     t = s[2..(size(s)-3)];
681  }
682  return(t);
683}
684
685static proc simfrac(poly up, poly down)
686{
687  // simplifies a fraction up/down
688  // into the form up/down = RT[1] + RT[2]/down
689  list LL = division(up,down);
690  list RT;
691  RT[1] =  LL[1][1,1]; // integer part
692  RT[2] = L[2][1]; // new numerator
693  return(RT);
694}
695
696proc findifs_example()
697"USAGE:  findifs_example();
698RETURN:  nothing (demo)
699PURPOSE:  demonstration of our approach and this library
700EXAMPLE: example findifs_example; shows examples
701"
702{
703
704  "* Equation: u_tt - A^2 u_xx -B^2 u_yy = 0; A,B are constants";
705  "* we employ three central differences";
706  "* the vector we act on is (u_xx, u_yy, u_tt, u)^T";
707  "* Set up the ring: ";
708  "ring r = (0,A,B,dt,dx,dy),(Tx,Ty,Tt),(c,dp);";
709  ring r = (0,A,B,dt,dx,dy),(Tx,Ty,Tt),(c,dp);
710  "* Set up the matrix with equation and approximations: ";
711  "matrix M[4][4] =";
712    "    // direct equation:";
713  "    -A^2, -B^2, 1, 0,";
714  "    // central difference u_tt";
715  "    0, 0,  -dt^2*Tt, (Tt-1)^2,";
716  "    // central difference u_xx";
717  "    -dx^2*Tx, 0, 0, (Tx-1)^2,";
718  "    // central difference u_yy";
719  "    0, -dy^2*Ty, 0, (Ty-1)^2;";
720  matrix M[4][4] =
721    // direct equation:
722    -A^2, -B^2, 1, 0,
723    // central difference u_tt
724    0, 0,  -dt^2*Tt, (Tt-1)^2,
725    // central difference u_xx
726    -dx^2*Tx, 0, 0, (Tx-1)^2,
727    // central difference u_yy
728    0, -dy^2*Ty, 0, (Ty-1)^2;
729//=========================================
730// CHECK THE CORRECTNESS OF EQUATIONS AS INPUT:
731ring rX = (0,A,B,dt,dx,dy,Tx,Ty,Tt),(Uxx, Uyy,Utt, U),(c,Dp);
732matrix M = imap(r,M);
733vector X = [Uxx, Uyy, Utt, U];
734 "* Print the differential form of equations: ";
735print(M*X);
736// END CHECK
737//=========================================
738 setring r;
739 "* Perform the elimination of module components:";
740 " module R = transpose(M);";
741 " module S = std(R);";
742 module R = transpose(M);
743 module S = std(R);
744 " * See the result of Groebner bases: generators are columns";
745 " print(S);";
746 print(S);
747 " * So, only the first column has its nonzero element in the last component";
748 " * Hence, this polynomial is the scheme";
749 " poly   p = S[4,1];" ;
750 poly   p = S[4,1]; // by elimination of module components
751 " print(p); ";
752 print(p);
753 list L; L[1]=A;L[2] = B;
754 ideal I = decoef(p,dt); // make splitting w.r.t. the appearance of dt
755 "* Create the nodal of the scheme in TeX format: ";
756 " ideal I = decoef(p,dt);";
757 " difpoly2tex(I,L);";
758 difpoly2tex(I,L); // the nodal form of the scheme in TeX
759 "* Preparations for the semi-factorized form: ";
760 poly pi1 = subst(I[2],B,0);
761 poly pi2 = I[2] - pi1;
762 " poly pi1 = subst(I[2],B,0);";
763 " poly pi2 = I[2] - pi1;";
764 "* Show the semi-factorized form of the scheme: 1st summand";
765 " factorize(I[1]); ";
766 factorize(I[1]); // semi-factorized form of the scheme: 1st summand
767 "* Show the semi-factorized form of the scheme: 2nd summand";
768 " factorize(pi1);";
769 factorize(pi1); // semi-factorized form of the scheme: 2nd summand
770 "* Show the semi-factorized form of the scheme: 3rd summand";
771 " factorize(pi1);";
772 factorize(pi2); // semi-factorized form of the scheme: 3rd summand
773}
774example
775{
776  "EXAMPLE:"; echo=1;
777 findifs_example();
778}
Note: See TracBrowser for help on using the repository browser.