# source:git/Singular/LIB/lejeune.lib@e7ff6b

spielwiese
Last change on this file since e7ff6b was e7ff6b, checked in by Nadine Cremer <cremer@…>, 18 years ago
*cremer: formaldiff seems to work,now start putting together... ;-) git-svn-id: file:///usr/local/Singular/svn/trunk@8384 2c84dea3-7e68-4137-9b89-c4e89433aadc
• Property mode set to 100644
File size: 3.7 KB
Line
1
2//-*- mode:C++;-*-
3// \$Id: lejeune.lib,v 1.4 2005-06-23 13:40:16 cremer Exp \$
4
5info="
6LIBRARY: lejeune1.4.lib  Arc space computations
8[SEE ALSO: <comma-separated words of cross references>]
9[KEYWORDS: <semicolon-separated phrases of index keys>]
10PROCEDURES:
11    variables(k,i);      creates k*i new var. t,a(1),..,a(i),..,x(1),..,x(i)
12    a_z(k);              returns kth letter of the alphabet
13    tpolys(k,i);         creates polyn. a(1)*t+..+a(n)*t^n
14    ringchange(i);       changes the ring to the one needed in ith step
15    plugin_coeffs(i,f)   plugs tpolys into f, up to power i
16    maxidealstep(i,N);   returns ideal needed for contraction in ith step;
17                         N is number of variables of input f
18    formaldiff(f,k);     computes the formal derivatives D_I with |I|<k
19  ";
20
21
22LIB "ring.lib";
23LIB "general.lib";
24
25
26
27
28proc formaldiff (poly f,int i,int a,int k)
29{
30  int s,t,v;                          // loop variables
31  int u;
32  def R=plugin_coeffs(f,i);           // plugs the power series in...
33  setring R;                          // changes the ring
34  def Coe=result;
35  matrix coe=Coe;                     // gives the t-coeff. after plugging in
36  poly fkv;                           // need this stuff for the following
37  ideal m;                            // loops...
38  ideal m1,m2,J,resultdiff;
39  for(v=0;v<=k-1;v++)                 // consider the different coeff.
40   {
41     fkv=coe[a+v,1];
42     m=fkv;
43     J=fkv;                           // will save the result in this step
44     for(s=1;s<k;s++)
45       {
46         m1=maxidealstep(i,startvar); // computes the corresponding ideal
47         m1=m1^s;
48         u=size(m1);
49         for(t=1;t<=u;t++)
50          {
51            m2=contract(m1[t],m);     // actual differentiation
52            J=J,m2;
53          }
54       }
55     resultdiff=resultdiff,J;
56   }
57  resultdiff;~
58  export(resultdiff);
59  return(R);
60
61}
62
63
64
65proc plugin_coeffs (poly f,int i)
66{
67  def r=basering;
68  def R=ringchange(i);
69  setring R;
70  ideal I=tpolys(i,startvar);
71  poly g=imap(r,f);
72  export(g);
73  map h=r,I;
74  ideal J=h(f);
75  export(h);
76  matrix result=coeffs(J[1],t);
77  export result;
78  return(R);
79}
80
81
82
83proc ringchange (int i)
84{
85  int startvar=nvars(basering);
86  export(startvar);
87  string str=variables(startvar,i);
88  def R=changevar(""+varstr(r)+",t,"+variables(startvar,i)+"");
89  return(R);
90}
91
92
93
94proc variables (int k,int i)
95{
96  list l;
97  int s,u;                              // loop variables
98  string str;
99  for (u=1;u<=k;u++)
100   {
101     for (s=1;s<=i;s++)
102     {
103       str=""+a_z(u)+"("+string(s)+")";
104       l[(u-1)*i+s]=str;
105     }
106   }
107  //l=insert(l,"t");
108  string str1=string(l);
109  return(str1);
110}
111
112
113proc a_z (int n)                      // returns nth letter of the alphabet
114{
115  if((n<1)||(n>26))                     // input admissible?
116   {
117     "n must range between 1 and 26!";
118      return(0);
119   }
120  string s="ring r=0,("+A_Z("a",n)+"),ds;";
121  execute(s);
122  return (string(var(n)));
123}
124
125
126
127proc tpolys (int i,int k)             // constructs polynomials a(1)*t+...
128{                                     // has to be called from tpolys
129  int s,t;                             // loop variables
130  int v;
131  poly sum;
132  ideal I;
133  for(t=1;t<=k;t++)
134   {
135     v=(t-1)*i;
136     for(s=1;s<=i;s++)
137      {
138        sum=sum+var(1+k+v+s)*var(k+1)^s;    // clumsy: working with "var(1)",
139      }                                     // depends on form of basering
140     I[t]=sum;
141     sum=0;
142   }
143  return(I);
144}
145
146
147
148
149proc maxidealstep (int i,int N)
150{
151  ideal I=var(N+1+i);
152  int j;
153  for(j=2;j<=N;j++)
154   {
155     I=I,var(N+1+j*i);
156   }
157return(I);
158}
Note: See TracBrowser for help on using the repository browser.