1 | //// |
---|
2 | version="version alexpoly.lib 4.0.0.0 Jun_2013 "; // $Id$ |
---|
3 | category="Singularities"; |
---|
4 | info=" |
---|
5 | LIBRARY: alexpoly.lib Resolution Graph and Alexander Polynomial |
---|
6 | AUTHOR: Fernando Hernando Carrillo, hernando@agt.uva.es |
---|
7 | Thomas Keilen, keilen@mathematik.uni-kl.de |
---|
8 | |
---|
9 | OVERVIEW: |
---|
10 | A library for computing the resolution graph of a plane curve singularity f, |
---|
11 | the total multiplicities of the total transforms of the branches of f along |
---|
12 | the exceptional divisors of a minimal good resolution of f, the Alexander |
---|
13 | polynomial of f, and the zeta function of its monodromy operator. |
---|
14 | |
---|
15 | PROCEDURES: |
---|
16 | resolutiongraph(f); resolution graph f |
---|
17 | totalmultiplicities(f); resolution graph, total multiplicities and strict multiplicities of f |
---|
18 | alexanderpolynomial(f); Alexander polynomial of f |
---|
19 | semigroup(f); calculates generators for the semigroup of f |
---|
20 | proximitymatrix(f); calculates the proximity matrix of f |
---|
21 | multseq2charexp(v); converts multiplicity sequence to characteristic exponents |
---|
22 | charexp2multseq(v); converts characteristic exponents to multiplicity sequence |
---|
23 | charexp2generators(v); converts characteristic exponents to generators of the semigroup |
---|
24 | charexp2inter(c,e); converts contact matrix and charact. exp. to intersection matrix |
---|
25 | charexp2conductor(v); converts characteristic exponents to conductor |
---|
26 | charexp2poly(v,a); calculates a polynomial f with characteristic exponents v |
---|
27 | tau_es2(f); equisingular Tjurina number of f |
---|
28 | |
---|
29 | KEYWORDS: Hamburger-Noether expansion; Puiseux expansion; curve singularities; |
---|
30 | topological invariants; Alexander polynomial; resolution graph; |
---|
31 | total multiplicities; equisingular Tjurina number |
---|
32 | "; |
---|
33 | |
---|
34 | /////////////////////////////////////////////////////////////////////////////////////////// |
---|
35 | LIB "hnoether.lib"; |
---|
36 | /////////////////////////////////////////////////////////////////////////////////////////// |
---|
37 | |
---|
38 | proc resolutiongraph (def INPUT,list #) |
---|
39 | "USAGE: resolutiongraph(INPUT); INPUT poly or list |
---|
40 | ASSUME: INPUT is either a REDUCED bivariate polynomial defining a plane curve singularity, |
---|
41 | or the output of @code{hnexpansion(f[,\"ess\"])}, or the list @code{hne} in |
---|
42 | the ring created by @code{hnexpansion(f[,\"ess\"])}, or the output of |
---|
43 | @code{develop(f)} resp. of @code{extdevelop(f,n)}, or a list containing |
---|
44 | the contact matrix and a list of integer vectors with the characteristic exponents |
---|
45 | of the branches of a plane curve singularity, or an integer vector containing the |
---|
46 | characteristic exponents of an irreducible plane curve singularity. |
---|
47 | RETURN: intmat, the incidence matrix of the resolution graph of the plane curve |
---|
48 | defined by INPUT, where the entries on the diagonal are the weights of the |
---|
49 | vertices of the graph and a negative entry corresponds to the strict transform |
---|
50 | of a branch of the curve. |
---|
51 | NOTE: In case the Hamburger-Noether expansion of the curve f is needed |
---|
52 | for other purposes as well it is better to calculate this first |
---|
53 | with the aid of @code{hnexpansion} and use it as input instead of |
---|
54 | the polynomial itself. |
---|
55 | @* |
---|
56 | If you are not sure whether the INPUT polynomial is reduced or not, use |
---|
57 | @code{squarefree(INPUT)} as input instead. |
---|
58 | SEE ALSO: develop, hnexpansion, totalmultiplicities, alexanderpolynomial |
---|
59 | EXAMPLE: example resolutiongraph; shows an example |
---|
60 | " |
---|
61 | { |
---|
62 | return(totalmultiplicities(INPUT,#)[1]); |
---|
63 | } |
---|
64 | example |
---|
65 | { |
---|
66 | "EXAMPLE:"; |
---|
67 | echo=2; |
---|
68 | ring r=0,(x,y),ls; |
---|
69 | poly f1=(y2-x3)^2-4x5y-x7; |
---|
70 | poly f2=y2-x3; |
---|
71 | poly f3=y3-x2; |
---|
72 | resolutiongraph(f1*f2*f3); |
---|
73 | } |
---|
74 | |
---|
75 | proc totalmultiplicities (def INPUT, list #) |
---|
76 | "USAGE: totalmultiplicities(INPUT); INPUT poly or list |
---|
77 | ASSUME: INPUT is either a REDUCED bivariate polynomial defining a plane curve singularity, |
---|
78 | or the output of @code{hnexpansion(f[,\"ess\"])}, or the list @code{hne} in |
---|
79 | the ring created by @code{hnexpansion(f[,\"ess\"])}, or the output of |
---|
80 | @code{develop(f)} resp. of @code{extdevelop(f,n)}, or a list containing |
---|
81 | the contact matrix and a list of integer vectors with the characteristic exponents |
---|
82 | of the branches of a plane curve singularity, or an integer vector containing the |
---|
83 | characteristic exponents of an irreducible plane curve singularity. |
---|
84 | RETURN: list @code{L} of three integer matrices. @code{L[1]} is the incidence matrix of |
---|
85 | the resolution graph of the plane curve defined by INPUT, where the entries on the |
---|
86 | diagonal are the weights of the vertices of the graph and a negative entry corresponds |
---|
87 | to the strict transform of a branch of the curve. @code{L[2]} is an integer matrix, |
---|
88 | which has for each vertex in the graph a row and for each branch of the curve a column. |
---|
89 | The entry in position [i,j] contains the total multiplicity of the j-th branch (i.e. the |
---|
90 | branch with weight -j in @code{L[1]}) along the exceptional divisor corresponding |
---|
91 | to the i-th row in @code{L[1]}. In particular, the i-th row contains |
---|
92 | the total multiplicities of the branches of the plane curve (defined by INPUT) along |
---|
93 | the exceptional divisor which corresponds to the i-th row in the incidence matrix |
---|
94 | @code{L[1]}. @code{L[3]} is an integer matrix which contains the (strict) multiplicities |
---|
95 | of the branches of the curve along the exceptional divisors in the same way as @code{L[2]} |
---|
96 | contains the total multiplicities. |
---|
97 | NOTE: The total multiplicty of a branch along an exceptional divisor is the multiplicity |
---|
98 | with which this exceptional divisor occurs in the total transform of this branch |
---|
99 | under the resolution corresponding to the resolution graph. |
---|
100 | @* |
---|
101 | In case the Hamburger-Noether expansion of the curve f is needed |
---|
102 | for other purposes as well it is better to calculate this first |
---|
103 | with the aid of @code{hnexpansion} and use it as input instead of |
---|
104 | the polynomial itself. |
---|
105 | @* |
---|
106 | If you are not sure whether the INPUT polynomial is reduced or not, use |
---|
107 | @code{squarefree(INPUT)} as input instead. |
---|
108 | SEE ALSO: develop, hnexpansion, alexanderpolynomial, resolutiongraph |
---|
109 | EXAMPLE: example totalmultiplicities; shows an example |
---|
110 | " |
---|
111 | { |
---|
112 | /////////////////////////////////////////////////////////////////////////////////////////////// |
---|
113 | /// The algorithm described in [JP] DeJong, Pfister, Local Analytic Geometry, Vieweg 2000, |
---|
114 | /// is used for the implementation -- the missing case is included by appropriate sorting. |
---|
115 | /////////////////////////////////////////////////////////////////////////////////////////////// |
---|
116 | /// If # is empty, then an additional sorting of the branches will be made, so that |
---|
117 | /// the branches can be split with split_graph in order to get the graphs of the curves |
---|
118 | /// separated on the first exceptional divisor. Otherwise split_graph should not be applied! |
---|
119 | /////////////////////////////////////////////////////////////////////////////////////////////// |
---|
120 | /// If #[1]="tau", then totalmultiplicities is called for the use in tau_es2 - thus it returns |
---|
121 | /// the prolonged resolutiongraphs of the branches, their multiplicity sequences and their |
---|
122 | /// contact matrix - sorted appropriately. |
---|
123 | /////////////////////////////////////////////////////////////////////////////////////////////// |
---|
124 | int i,j,s,e,ii; |
---|
125 | intmat helpmati,helpmatii; |
---|
126 | intvec helpvec; |
---|
127 | ///////////////////////////////////////////////////////////////////////////////// |
---|
128 | // Decide, which kind of input we have, and define the contact matrix |
---|
129 | ///////////////////////////////////////////////////////////////////////////////// |
---|
130 | // Input: a polynomial defining a plane curve singularity. |
---|
131 | ////////////////////////////////////////////////////////////////////////////// |
---|
132 | if (typeof(INPUT)=="poly") |
---|
133 | { |
---|
134 | /* |
---|
135 | poly f=squarefree(INPUT); |
---|
136 | if ( deg(f)!=deg(INPUT) ) |
---|
137 | { |
---|
138 | dbprint(printlevel-voice+4,"// input polynomial was not reduced"); |
---|
139 | dbprint(printlevel-voice+4,"// we continue with its reduction"); |
---|
140 | } |
---|
141 | list I@N@V=invariants(f); |
---|
142 | */ |
---|
143 | list I@N@V=invariants(INPUT); |
---|
144 | } |
---|
145 | else |
---|
146 | { |
---|
147 | /////////////////////////////////////////////////////////////////////////////////// |
---|
148 | // Input: the vector of characteristic exponents of an irreducible plane curve |
---|
149 | /////////////////////////////////////////////////////////////////////////////////// |
---|
150 | if (typeof(INPUT)=="intvec") |
---|
151 | { |
---|
152 | list charexp; |
---|
153 | charexp[1]=INPUT; |
---|
154 | intmat contact[1][1]=0; |
---|
155 | } |
---|
156 | else |
---|
157 | { |
---|
158 | if (typeof(INPUT)=="list") |
---|
159 | { |
---|
160 | ///////////////////////////////////////////////////////////////////////////////// |
---|
161 | // Input: intersection-matrix and characteristic exponents. |
---|
162 | ////////////////////////////////////////////////////////////////////////////// |
---|
163 | if (typeof(INPUT[1])=="intmat") |
---|
164 | { |
---|
165 | intmat contact=INPUT[1]; |
---|
166 | list charexp=INPUT[2]; |
---|
167 | } |
---|
168 | else |
---|
169 | { |
---|
170 | ///////////////////////////////////////////////////////////////////////////////// |
---|
171 | // Input: output of hnexpansion or hne in the ring created by hnexpansion |
---|
172 | ////////////////////////////////////////////////////////////////////////////// |
---|
173 | if ((typeof(INPUT[1])=="ring") or (typeof(INPUT[1])=="list")) |
---|
174 | { |
---|
175 | list I@N@V=invariants(INPUT); |
---|
176 | } |
---|
177 | else |
---|
178 | { |
---|
179 | //////////////////////////////////////////////////////////////////////////////////// |
---|
180 | // Input: output of reddevelop or extdevelop -- irreducible plane curve singularity |
---|
181 | //////////////////////////////////////////////////////////////////////////////////// |
---|
182 | if (typeof(INPUT[1])=="matrix") |
---|
183 | { |
---|
184 | list charexp=invariants(INPUT)[1]; |
---|
185 | intmat contact[1][1]=0; |
---|
186 | } |
---|
187 | else |
---|
188 | { |
---|
189 | ERROR("The input is invalid!"); |
---|
190 | } |
---|
191 | } |
---|
192 | } |
---|
193 | } |
---|
194 | else |
---|
195 | { |
---|
196 | ERROR("The input is invalid!"); |
---|
197 | } |
---|
198 | } |
---|
199 | } |
---|
200 | /////////////////////////////////////////////////////////////////////////////////////////////////// |
---|
201 | // If the input was a polynomial or a HN-Expansion, then calculate the contact matrix and char.exponents. |
---|
202 | /////////////////////////////////////////////////////////////////////////////////////////////////// |
---|
203 | if (defined(I@N@V)) |
---|
204 | { |
---|
205 | intmat contact=I@N@V[size(I@N@V)][1]; // contact numbers |
---|
206 | list charexp; // characteristic exponents |
---|
207 | for (i=1;i<=size(I@N@V)-1;i++) |
---|
208 | { |
---|
209 | charexp[i]=I@N@V[i][1]; |
---|
210 | } |
---|
211 | } |
---|
212 | ////////////////////////////////////////////////////////////////////////////// |
---|
213 | // Find the maximal contact numbers of the branches |
---|
214 | ////////////////////////////////////////////////////////////////////////////// |
---|
215 | int r=ncols(contact); // number of branches |
---|
216 | intvec maxcontact; // maximal contactnumber of branch i with other branches |
---|
217 | for (i=1;i<=r;i++) |
---|
218 | { |
---|
219 | maxcontact[i]=max_in_intvec(intvec(contact[i,1..r])); |
---|
220 | } |
---|
221 | /////////////////////////////////////////////////////////////////////////////// |
---|
222 | // Define the graphs of the branches and prolong them if necessary. |
---|
223 | /////////////////////////////////////////////////////////////////////////////// |
---|
224 | intmat gr,gr_red,grp; // a subgraph, a reduced subgraph, a prolonged subgraph |
---|
225 | int omega; // point of highest weight in subgraph |
---|
226 | list graphs; // contains the subgraphs of the C_i |
---|
227 | list totmult; // contains the total multiplicities of the subgraphs |
---|
228 | list multipl; // contains the multiplicities of the subgraphs |
---|
229 | list gr_tm; // takes a single subgraph and tot.mult. |
---|
230 | intvec tm,mt; // total multiplicities and multiplicities |
---|
231 | for (i=1;i<=r;i++) |
---|
232 | { |
---|
233 | gr_tm=irred_resgraph_totmult(charexp[i]); // graph, total multipl. and multipl. of the ith branch |
---|
234 | gr=gr_tm[1]; // the graph of the ith branch |
---|
235 | tm=gr_tm[2]; // the total multiplicities of the ith branch |
---|
236 | mt=gr_tm[3]; // the multiplicities of the ith branch |
---|
237 | omega=nrows(gr)-1; // maximal weight in the graph of the ith branch |
---|
238 | // If the maximal contact of the ith branch with some other branch is larger |
---|
239 | // than the maximal weight omega, then the graph has to be polonged. |
---|
240 | if (omega<maxcontact[i]) |
---|
241 | { |
---|
242 | grp=intmat(intvec(0),maxcontact[i]+1,maxcontact[i]+1); |
---|
243 | // save the graph without the point of the strict transform |
---|
244 | if (omega>=1) // otherwise the graph consists only of the strict transform |
---|
245 | { |
---|
246 | grp[1..omega,1..omega]=gr[1..omega,1..omega]; |
---|
247 | } |
---|
248 | // add the points of multiplicity 1 up to weight maxcontact[i] and the strict transform |
---|
249 | // and connect them |
---|
250 | for (j=omega+1;j<=maxcontact[i]+1;j++) |
---|
251 | { |
---|
252 | // adding the vertex to the graph and adding the total multiplicities |
---|
253 | if (j<=maxcontact[i]) |
---|
254 | { |
---|
255 | grp[j,j]=j; |
---|
256 | if (j>1) |
---|
257 | { |
---|
258 | tm[j]=tm[j-1]+1; |
---|
259 | mt[j]=1; |
---|
260 | } |
---|
261 | else // then there is no previous point in the graph |
---|
262 | { |
---|
263 | tm[j]=1; |
---|
264 | mt[j]=1; |
---|
265 | } |
---|
266 | } |
---|
267 | // connecting the vertex with the previous part of the graph |
---|
268 | if (j>1) // otherwise there is nothing to connect to |
---|
269 | { |
---|
270 | grp[j-1,j]=1; |
---|
271 | grp[j,j-1]=1; |
---|
272 | } |
---|
273 | } |
---|
274 | gr=grp; // replace the subgraph by the prolonged subgraph |
---|
275 | } |
---|
276 | gr[nrows(gr),ncols(gr)]=-i; // give the strict transform in the ith branch weight -i |
---|
277 | graphs[i]=gr; |
---|
278 | totmult[i]=tm; |
---|
279 | multipl[i]=mt; |
---|
280 | } |
---|
281 | /////////////////////////////////////////////////////////////////////////////// |
---|
282 | // Check, if the procedure is called from inside tau_es2. |
---|
283 | int check_tau; |
---|
284 | if (size(#)==1) |
---|
285 | { |
---|
286 | if (#[1]=="tau") |
---|
287 | { |
---|
288 | check_tau=1; |
---|
289 | } |
---|
290 | } |
---|
291 | ///////////////////////////////////////////////////////////////////////////////// |
---|
292 | // The procedure tau_es2 expects the branches to be ordered according to their |
---|
293 | // contact during the resolution process. |
---|
294 | ///////////////////////////////////////////////////////////////////////////////// |
---|
295 | if ((size(#)==0) or (check_tau==1)) |
---|
296 | { |
---|
297 | list SORT; |
---|
298 | for (i=1;i<=ncols(contact)-2;i++) |
---|
299 | { |
---|
300 | SORT=sort_branches(contact,graphs,totmult,multipl,i,ncols(contact)); |
---|
301 | contact=SORT[1]; |
---|
302 | graphs=SORT[2]; |
---|
303 | totmult=SORT[3]; |
---|
304 | multipl=SORT[4]; |
---|
305 | } |
---|
306 | } |
---|
307 | /////////////////////////////////////////////////////////////////////////////////// |
---|
308 | /// If the procedure is called from inside tau_es2, the output will be the prolonged |
---|
309 | /// resolutiongraphs of the branches, their multiplicity sequences and their contact matrix. |
---|
310 | if (check_tau==1) |
---|
311 | { |
---|
312 | return(list(graphs,multipl,contact)); |
---|
313 | } |
---|
314 | ///////////////////////////////////////////////////////////////////////////////////// |
---|
315 | // Sort the branches in such a way, that in the blowing up procedure whenever to branches |
---|
316 | // C_i and C_j have contact k (i.e. after k steps they are separated) and one of the following |
---|
317 | // situations occurs: |
---|
318 | // A) some exceptional divisor E_l (l<k) still intersects C_i after k steps of blowing up |
---|
319 | // and no such E_l intersects C_j |
---|
320 | // OR |
---|
321 | // B) some exceptional divisor E_l (l<k-1) still intersects C_i after k steps of blowing up |
---|
322 | // and another one (necessarily E_k-1) intersects C_j, |
---|
323 | // THEN: i<j! |
---|
324 | ///////////////////////////////////////////////////////////////////////////////////////// |
---|
325 | // This means in the algorithm for glueing graphs together described in [JP] p.217-19 |
---|
326 | // the Case 3 never occurs -- and Case 1 only occurs, if it is unavoidable, that is when |
---|
327 | // C_1 there meets an E_l with l<k. |
---|
328 | // (Otherwise we say that (C_i,C_j) has "bad contact".) |
---|
329 | //////////////////////////////////////////////////////////////////////////////////////// |
---|
330 | // We can detect this by looking at the graphs of C_i and C_j. If E_l stays together with |
---|
331 | // C_i and not with C_j in the k-th step of blowing up, then in the graph of C_i |
---|
332 | // k and l are NOT connected, while in the graph of C_j they are! |
---|
333 | ////////////////////////////////////////////////////////////////////////////////////////// |
---|
334 | // Note also, that for a fixed pair (i,j) this situation can only occur for ONE l (in Case A: l<k |
---|
335 | // and possibly l=k-1; in Case B: l<k-1). |
---|
336 | ////////////////////////////////////////////////////////////////////////////////////////// |
---|
337 | // Our algorithm now works as follows. Start with i=1 and j=2. |
---|
338 | // While (i<r=number of branches) do as follows: Check C_i with C_j. |
---|
339 | // If a bad contact (C_i,C_j) occurs, then MOVE C_j to the position of C_i (and move all C_t |
---|
340 | // with t>=i, t!=j one position further). Else do nothing particular. |
---|
341 | // Then, if j<r set j=j+1. Else set j=i+2 and i=i+1. End of While. |
---|
342 | // To see that this works we note that if (C_i,C_j) has bad contact for some j>i, then |
---|
343 | // (C_j,C_t) does NOT have bad contact for t=i+1,...,j-1. For this we consider the 3 cases: |
---|
344 | // 1) contact(C_t,C_j)=contact(C_i,C_j)=k: E_l stays with C_j in k-th step and C_t and C_j |
---|
345 | // separate there, so (C_t,C_j) has bad contact |
---|
346 | // and hence (C_j,C_t) does NOT have bad contact |
---|
347 | // 2) contact(C_t,C_j)>contact(C_i,C_j)=k: CANNOT HAPPEN - C_i had no bad contact with C_t, |
---|
348 | // hence E_l did not stay with C_t in the k-th step; |
---|
349 | // however, E_l stays with C_j there, so C_t and C_j |
---|
350 | // cannot have a contact higher than k |
---|
351 | // 3) contact(C_t,C_j)<contact(C_i,C_j)=k: when C_t and C_j separate, C_j and C_i are still |
---|
352 | // together; and since (C_i,C_t) did not have bad |
---|
353 | // contact also (C_j,C_t) cannot have bad contact. |
---|
354 | // This ensures that when we do the moving of the graphs, we do NOT create any new bad contacts. |
---|
355 | ////////////////////////////////////////////////////////////////////////////// |
---|
356 | i=1; |
---|
357 | j=2; |
---|
358 | intvec connections_i,connections_j; |
---|
359 | int n_i,n_j; |
---|
360 | int bad_contact; |
---|
361 | while (i<r) |
---|
362 | { |
---|
363 | // Test graphs[i] with graphs[j] for bad contact. |
---|
364 | connections_i=find_lower_connection_points(graphs[i],contact[i,j]); |
---|
365 | if (connections_i[1]==0) // no connection point of lower weight there |
---|
366 | { |
---|
367 | n_i=0; |
---|
368 | } |
---|
369 | else |
---|
370 | { |
---|
371 | n_i=size(connections_i); |
---|
372 | } |
---|
373 | connections_j=find_lower_connection_points(graphs[j],contact[i,j]); |
---|
374 | if (connections_j[1]==0) // no connection point of lower weight there |
---|
375 | { |
---|
376 | n_j=0; |
---|
377 | } |
---|
378 | else |
---|
379 | { |
---|
380 | n_j=size(connections_j); |
---|
381 | } |
---|
382 | bad_contact=0; |
---|
383 | ii=1; |
---|
384 | // The points of weight k=contact(C_i,C_j) in C_i and C_j are connectd to n_i respectively |
---|
385 | // n_j points of weight less than k, where n_i and n_j might take values among 0, 1 and 2. |
---|
386 | // n_j=2: Then E_k is intersected by E_k-1 and E_l (l<k-1) in the k-th step, however |
---|
387 | // C_j does not stay with anyone of those, so (C_i,C_j) does not have bad contact. |
---|
388 | // n_i=0: Then C_i stays with E_k-1 and there is no E_l with l<k-1 in the k-th step, |
---|
389 | // hence (C_i,C_j) does not have bad contact. |
---|
390 | // n_i=1, n_j=0: Analogously, then C_j stays with E_k-1 and there is no E_l (l<k-1). This means |
---|
391 | // (C_i,C_j) has bad contact. |
---|
392 | // n_i=1, n_j=1: EITHER there is no E_l (l<k-1) and C_i and C_j both separate from E_k-1 |
---|
393 | // in the k-th step of blowing up (hence (C_i,C_j) does not have bad contact). |
---|
394 | // OR there is an E_l (l<k-1) with which one stays and the other one stays |
---|
395 | // with E_k-1 (bad contact, if E_l stays with C_j; good contact otherwise). |
---|
396 | // n_i=2, E_k is intersected by E_k-1 and E_l in the k-th step and C_i does not |
---|
397 | // stay together with any of those. |
---|
398 | // n_j=0: This CANNOT occur, since then C_j would have to stay together with E_l |
---|
399 | // and E_k-1, which is impossible. |
---|
400 | // n_j=1: Then C_j stays together with either E_l or with E_k-1, however both cases |
---|
401 | // imply that (C_i,C_j) has bad contact |
---|
402 | // ALTERNATIVE CONSIDERATION FOR DISTINGUISHING THE CASES: |
---|
403 | // n_i<n_j : C_i stays with "more" E's than C_j does (only one is possible of course!), |
---|
404 | // hence (C_i,C_j) does not have bad contact. |
---|
405 | // n_i>n_j : C_i stays with "less" E's than C_j does (i.e. with none) hence (C_i,C_j) has bad contact. |
---|
406 | // n_i=n_j=0 : Not possible, since then both would stay with E_k-1. |
---|
407 | // n_i=n_j=2 : Both separate from E_k-1 as well as from some E_l (l<k) - hence NO bad contact. |
---|
408 | // n_i=n_j=1 : One stays with E_k-1 and one with E_l (l<k). Bad contact, if C_j stays with E_l, |
---|
409 | // i.e. if in the graph of C_i the point k is connected to l, and no bad contact otherwise. |
---|
410 | if ((n_i>n_j) or ((n_i==1) and (n_j==1) and (connections_i[1]<contact[i,j]-1))) |
---|
411 | { |
---|
412 | bad_contact=1; |
---|
413 | // Move the graph (etc.) of C_j into the position i. |
---|
414 | graphs=insert(graphs,graphs[j],i-1); |
---|
415 | graphs=delete(graphs,j+1); |
---|
416 | totmult=insert(totmult,totmult[j],i-1); |
---|
417 | totmult=delete(totmult,j+1); |
---|
418 | multipl=insert(multipl,multipl[j],i-1); |
---|
419 | multipl=delete(multipl,j+1); |
---|
420 | contact=move_row_col(contact,i,j); |
---|
421 | } |
---|
422 | if (j<r) // There are still some C_j against which C_i has to be tested. |
---|
423 | { |
---|
424 | j++; |
---|
425 | } |
---|
426 | else // Now C_i has been tested against all C_j, and we may continue with C_i+1. |
---|
427 | { |
---|
428 | j=i+2; |
---|
429 | i=i+1; |
---|
430 | } |
---|
431 | } |
---|
432 | ///////////////////////////////////////////////////////////////////////////////////// |
---|
433 | // Glue the graphs together. |
---|
434 | ////////////////////////////////////////////////////////////////////////////// |
---|
435 | /////////////////////////////////////////////////////////////////////////////////// |
---|
436 | intmat rgraph=graphs[1]; // keeps the resolution graph |
---|
437 | intmat rtm=intmat(totmult[1],nrows(rgraph),1); // keeps the tot.mult. at the vertices of the graph as rows |
---|
438 | intmat rmt=intmat(multipl[1],nrows(rgraph),1); // keeps the mult. at the vertices of the graph as rows |
---|
439 | intvec stricttransforms=0,ncols(rgraph); // keeps the position of the ith strict |
---|
440 | // transform in rgraph at position i+1 !!! |
---|
441 | intvec k,kp,p,q,o; // highest contact numbers, num. of branch with hgt contact, separation points |
---|
442 | int maxc,maxcp; |
---|
443 | for (i=2;i<=r;i++) |
---|
444 | { |
---|
445 | ////////////////////////////////////////////////////////////////////////////////// |
---|
446 | // Find j<i minimal s.t. contact[i,j] is maximal. |
---|
447 | maxcp=i; |
---|
448 | for (j=1;j<i;j++) |
---|
449 | { |
---|
450 | if (contact[i,j]>contact[i,maxcp]){maxcp=j;} |
---|
451 | } |
---|
452 | kp[i]=maxcp; // the j such that C_i has its maximal contact to C_j |
---|
453 | k[i]=contact[i,maxcp]; // the maximal contact of C_i with C_1,...,C_i-1 |
---|
454 | /////////////////////////////////////////////////////////////////////////////////// |
---|
455 | // Find in the graph of C_1,...,C_i-1 the points p of wgt k and q of wgt k-1 |
---|
456 | // connected to C_maxcp. |
---|
457 | // Since non of C_j for j<maxcp has contact k with C_i, the point p lies in |
---|
458 | // the remaining part of the graph of C_maxcp. |
---|
459 | s=rgraph[stricttransforms[maxcp]+1,stricttransforms[maxcp]+1]; |
---|
460 | p[i]=stricttransforms[maxcp]+1+k[i]-s; // pt. to which reduced subgraph of C_i is glued |
---|
461 | // If s<k[i], then q also lies in this part, otherwise it lies in the remaining part |
---|
462 | // of the graph of the C_j to which C_maxcp is connected, i.e. j=kp[maxcp], since |
---|
463 | // the contact of C_i and of C_maxcp to this C_j is strictly less than k. |
---|
464 | // If s=k[i]=1, then p=1 and there is no q! We may thus set q=0. |
---|
465 | if ((s<k[i]) or ((s==k[i]) and (s==1))) // i.e. q is on the same subgraph as p, or q does not exist |
---|
466 | { |
---|
467 | q[i]=p[i]-1; |
---|
468 | } |
---|
469 | else // i.e. q is on the subgraph to which the subgraph of p has been glued |
---|
470 | { |
---|
471 | s=rgraph[stricttransforms[kp[maxcp]]+1,stricttransforms[kp[maxcp]]+1]; |
---|
472 | q[i]=stricttransforms[kp[maxcp]]+k[i]-s; |
---|
473 | } |
---|
474 | ////////////////////////////////////////////////////////////////////////////////////// |
---|
475 | // Reduce the graph of C_i and add it to the graph of C_1,...,C_i-1. |
---|
476 | gr=graphs[i]; |
---|
477 | s=nrows(rgraph); |
---|
478 | // Delete in the graph of C_i everything of weight <=k. |
---|
479 | gr_red=intmat(intvec(gr[k[i]+1..nrows(gr),k[i]+1..ncols(gr)]),nrows(gr)-k[i],ncols(gr)-k[i]); |
---|
480 | // Add this part to the graph of C_1,...,C_i-1. |
---|
481 | rgraph=addmat(rgraph,gr_red); |
---|
482 | ///////////////////////////////////////////////////////////////////////////////////// |
---|
483 | // Connect the two parts of the graph. |
---|
484 | ///////////////////////////////////////////////////////////////////////////////////// |
---|
485 | // Connect the points connected to the point of wgt k in the graph of C_i to p[i]. |
---|
486 | for (j=k[i]+1;j<=ncols(gr);j++) |
---|
487 | { |
---|
488 | if(gr[k[i],j]==1) |
---|
489 | { |
---|
490 | rgraph[s+j-k[i],p[i]]=1; |
---|
491 | rgraph[p[i],s+j-k[i]]=1; |
---|
492 | } |
---|
493 | } |
---|
494 | // If pt. of wgt k is not connected to pt of wgt k-1 in graph of C_i, then points |
---|
495 | // connected to the one of wgt k-1 have to be connected to q[i]. |
---|
496 | if (k[i]>1) |
---|
497 | { |
---|
498 | if (gr[k[i],k[i]-1]!=1) |
---|
499 | { |
---|
500 | for (j=k[i]+1;j<=ncols(gr);j++) |
---|
501 | { |
---|
502 | if(gr[k[i]-1,j]==1) |
---|
503 | { |
---|
504 | rgraph[s+j-k[i],q[i]]=1; |
---|
505 | rgraph[q[i],s+j-k[i]]=1; |
---|
506 | } |
---|
507 | } |
---|
508 | // Cut the connection from p[i] to q[i]. |
---|
509 | rgraph[p[i],q[i]]=0; |
---|
510 | rgraph[q[i],p[i]]=0; |
---|
511 | } |
---|
512 | } |
---|
513 | stricttransforms[i+1]=ncols(rgraph); |
---|
514 | //////////////////////////////////////////////////////////////////////////////// |
---|
515 | // Adjust the total multiplicities |
---|
516 | //////////////////////////////////////////////////////////////////////////////// |
---|
517 | // Add the total multiplicities for the added subgraph to rtm |
---|
518 | tm=totmult[i]; |
---|
519 | mt=multipl[i]; |
---|
520 | if (k[i]<size(tm)) // if the reduced subgraph of C_i has more than one point |
---|
521 | { |
---|
522 | rtm=addmat(rtm,intmat(intvec(tm[k[i]+1..size(tm)]),nrows(gr_red),1)); |
---|
523 | rmt=addmat(rmt,intmat(intvec(mt[k[i]+1..size(mt)]),nrows(gr_red),1)); |
---|
524 | } |
---|
525 | else // if the reduced subgraph of C_i consists only of the strict transform |
---|
526 | { |
---|
527 | rtm=addmat(rtm,0); |
---|
528 | rmt=addmat(rmt,0); |
---|
529 | } |
---|
530 | // Adjust the total multiplicities at the places where the subgraph has been glued. |
---|
531 | e=k[i]; // the highest weight of a point that has not yet been assigned its tot. mult. |
---|
532 | while(e>=1) |
---|
533 | { |
---|
534 | s=stricttransforms[maxcp]+1; // Line in the graph of the start. pt. of the subgraph of C_maxcp. |
---|
535 | for (j=rgraph[s,s];j<=e;j++) // Adjust the multiplicities. |
---|
536 | { |
---|
537 | rtm[s+j-rgraph[s,s],i]=tm[j]; |
---|
538 | rmt[s+j-rgraph[s,s],i]=mt[j]; |
---|
539 | } |
---|
540 | maxcp=kp[maxcp]; // To which subgraph has the subgraph of C_maxcp been glued? |
---|
541 | e=rgraph[s,s]-1; // What is the highest weight of a pt. that has not yet been assigned its tot.mult.? |
---|
542 | } |
---|
543 | e=nrows(rtm); // Number of rows in the matrix of totalmultiplicities. |
---|
544 | // The total multiplicities of the C_s for s=1,...,i-1 along the exceptional divisors |
---|
545 | // which are introduced after the strict transform of C_s has separated (i.e. the entries |
---|
546 | // in rows stricttransform[i]+1,...,stricttransform[i+1]-1 in the s-th column of the matrix |
---|
547 | // of total multiplicities still have to be calculated. |
---|
548 | for (s=1;s<i;s++) |
---|
549 | { |
---|
550 | rtm[1..e,s]=adjust_tot_mult(intvec(rtm[1..e,i]),intvec(rtm[1..e,s]),intvec(rmt[1..e,i]),intvec(rmt[1..e,s]),p,q,stricttransforms,i); |
---|
551 | } |
---|
552 | // The total multiplicities of the C_i along the exceptional divisors |
---|
553 | // which are introduced for the sake of C_s, s=1,...,i-1, after the strict transform |
---|
554 | // of C_i has separated (i.e. the entries in rows stricttransform[s]+1,...,stricttransform[s+1]-1 |
---|
555 | // in the i-th column of the matrix of total multiplicities still have to be calculated. |
---|
556 | for (s=1;s<i;s++) |
---|
557 | { |
---|
558 | rtm[1..e,i]=adjust_tot_mult(intvec(rtm[1..e,s]),intvec(rtm[1..e,i]),intvec(rmt[1..e,s]),intvec(rmt[1..e,i]),p,q,stricttransforms,s); |
---|
559 | } |
---|
560 | } |
---|
561 | list result=rgraph,rtm,rmt; |
---|
562 | return(result); |
---|
563 | } |
---|
564 | example |
---|
565 | { |
---|
566 | "EXAMPLE:"; |
---|
567 | echo=2; |
---|
568 | ring r=0,(x,y),ls; |
---|
569 | poly f1=(y2-x3)^2-4x5y-x7; |
---|
570 | poly f2=y2-x3; |
---|
571 | poly f3=y3-x2; |
---|
572 | totalmultiplicities(f1*f2*f3); |
---|
573 | } |
---|
574 | |
---|
575 | |
---|
576 | |
---|
577 | proc alexanderpolynomial (def INPUT) |
---|
578 | "USAGE: alexanderpolynomial(INPUT); INPUT poly or list |
---|
579 | ASSUME: INPUT is either a REDUCED bivariate polynomial defining a plane curve singularity, |
---|
580 | or the output of @code{hnexpansion(f[,\"ess\"])}, or the list @code{hne} in |
---|
581 | the ring created by @code{hnexpansion(f[,\"ess\"])}, or the output of |
---|
582 | @code{develop(f)} resp. of @code{extdevelop(f,n)}, or a list containing |
---|
583 | the contact matrix and a list of integer vectors with the characteristic exponents |
---|
584 | of the branches of a plane curve singularity, or an integer vector containing the |
---|
585 | characteristic exponents of an irreducible plane curve singularity. |
---|
586 | CREATE: a ring with variables t, t(1), ..., t(r) (where r is the number of branches of |
---|
587 | the plane curve singularity f defined by INPUT) and ordering ls over the |
---|
588 | ground field of the basering. @* |
---|
589 | Moreover, the ring contains the Alexander polynomial of f in variables t(1), ..., t(r) |
---|
590 | (@code{alexpoly}), the zeta function of the monodromy operator of f in the variable t |
---|
591 | (@code{zeta_monodromy}), and a list containing the factors of the Alexander |
---|
592 | polynomial with multiplicities (@code{alexfactors}). |
---|
593 | RETURN: a list, say @code{ALEX}, where @code{ALEX[1]} is the created ring |
---|
594 | NOTE: to use the ring type: @code{def ALEXring=ALEX[i]; setring ALEXring;}. |
---|
595 | @* |
---|
596 | Alternatively you may use the procedure sethnering and type: sethnering(ALEX,\"ALEXring\"); |
---|
597 | @* |
---|
598 | To access the Alexander polynomial resp. the zeta function resp. the |
---|
599 | factors of the Alexander polynomial type: @code{alexpoly} resp. @code{zeta_monodromy} |
---|
600 | resp. @code{alexfactors}.@* |
---|
601 | In case the Hamburger-Noether expansion of the curve f is needed |
---|
602 | for other purposes as well it is better to calculate this first |
---|
603 | with the aid of @code{hnexpansion} and use it as input instead of |
---|
604 | the polynomial itself. |
---|
605 | @* |
---|
606 | If you are not sure whether the INPUT polynomial is reduced or not, use |
---|
607 | @code{squarefree(INPUT)} as input instead. |
---|
608 | SEE ALSO: resolutiongraph, totalmultiplicities |
---|
609 | EXAMPLE: example alexanderpolynomial; shows an example |
---|
610 | " |
---|
611 | { |
---|
612 | // Get the resolution graph and the total multiplicities. |
---|
613 | list gr_tm=totalmultiplicities(INPUT); |
---|
614 | intmat gr=gr_tm[1]; |
---|
615 | intmat tm=gr_tm[2]; |
---|
616 | int r=ncols(tm); |
---|
617 | int e=ncols(gr); |
---|
618 | // Define the Ring for the Alexander Polynomial and the Zeta Function of the Monodromy. |
---|
619 | execute("ring ALEXring=("+charstr(basering)+"),(t,t(1..r)),dp;"); |
---|
620 | poly hilfspoly=1; |
---|
621 | poly alexnumerator=1; // numerator of the Alexander polynomial |
---|
622 | poly alexdenominator=1; // denominator of the Alexander polynomial |
---|
623 | list alexfactors; // the factors of the Alexanderpolynomial with multiplicities |
---|
624 | list alexfactor; |
---|
625 | int chi=2; |
---|
626 | int i,j,k; |
---|
627 | int s=1; |
---|
628 | // Consider every vertex of the resolution graph. |
---|
629 | for (i=1;i<=e;i++) |
---|
630 | { |
---|
631 | if (gr[i,i]>0) // If it belongs to an exceptional curve. |
---|
632 | { |
---|
633 | for (j=1;j<=e;j++) // Calculate the Euler charateristik of the smooth locus of the exc. curve. |
---|
634 | { |
---|
635 | if ((gr[i,j]==1) and (i!=j)) |
---|
636 | { |
---|
637 | chi=chi-1; |
---|
638 | } |
---|
639 | } |
---|
640 | if (chi!=0) // If the Euler characteristik is not zero, then it gives a factor in the AP. |
---|
641 | { |
---|
642 | for (k=1;k<=r;k++) |
---|
643 | { |
---|
644 | hilfspoly=hilfspoly*t(k)^tm[i,k]; |
---|
645 | } |
---|
646 | hilfspoly=1-hilfspoly; |
---|
647 | if (chi<0) // ... either in the numerator ... |
---|
648 | { |
---|
649 | alexnumerator=alexnumerator * hilfspoly^-chi; |
---|
650 | } |
---|
651 | else // ... or in the denominator. |
---|
652 | { |
---|
653 | alexdenominator=alexdenominator * hilfspoly^chi; |
---|
654 | } |
---|
655 | alexfactor=hilfspoly,-chi; |
---|
656 | alexfactors[s]=alexfactor; |
---|
657 | s++; |
---|
658 | } |
---|
659 | chi=2; |
---|
660 | hilfspoly=1; |
---|
661 | } |
---|
662 | } |
---|
663 | // Calculate the Alexander Polynomial. |
---|
664 | if (ncols(tm)==1) // If we have only one branch, then the numerator has to be multiplied by 1-t. |
---|
665 | { |
---|
666 | alexnumerator=alexnumerator*(1-t(1)); |
---|
667 | alexfactor=1-t(1),1; |
---|
668 | alexfactors[size(alexfactors)+1]=alexfactor; |
---|
669 | } |
---|
670 | poly alexpoly=alexnumerator / alexdenominator; |
---|
671 | // Calculate the Zeta Function of the Monodromy Operator. |
---|
672 | poly zeta_monodromy=alexpoly; |
---|
673 | for (i=1;i<=r;i++) |
---|
674 | { |
---|
675 | zeta_monodromy=subst(zeta_monodromy,t(i),t); |
---|
676 | } |
---|
677 | export zeta_monodromy; |
---|
678 | export alexnumerator; |
---|
679 | export alexdenominator; |
---|
680 | export alexfactors; |
---|
681 | export alexpoly; |
---|
682 | list ALEX=ALEXring; |
---|
683 | return(ALEX); |
---|
684 | } |
---|
685 | example |
---|
686 | { |
---|
687 | "EXAMPLE:"; |
---|
688 | echo=2; |
---|
689 | ring r=0,(x,y),ls; |
---|
690 | poly f1=(y2-x3)^2-4x5y-x7; |
---|
691 | poly f2=y2-x3; |
---|
692 | poly f3=y3-x2; |
---|
693 | list ALEX=alexanderpolynomial(f1*f2*f3); |
---|
694 | def ALEXring=ALEX[1]; |
---|
695 | setring ALEXring; |
---|
696 | alexfactors; |
---|
697 | alexpoly; |
---|
698 | zeta_monodromy; |
---|
699 | } |
---|
700 | |
---|
701 | proc semigroup (def INPUT) |
---|
702 | "USAGE: semigroup(INPUT); INPUT poly or list |
---|
703 | ASSUME: INPUT is either a REDUCED bivariate polynomial defining a plane curve singularity, |
---|
704 | or the output of @code{hnexpansion(f[,\"ess\"])}, or the list @code{hne} in |
---|
705 | the ring created by @code{hnexpansion(f[,\"ess\"])}, or the output of |
---|
706 | @code{develop(f)} resp. of @code{extdevelop(f,n)}, or a list containing |
---|
707 | the contact matrix and a list of integer vectors with the characteristic exponents |
---|
708 | of the branches of a plane curve singularity, or an integer vector containing |
---|
709 | the characteristic exponents of an irreducible plane curve singularity. |
---|
710 | RETURN: a list with three entries. The first and the second are lists @code{v_1,...,v_s} |
---|
711 | and @code{w_1,...,w_r} respectively of integer vectors such that the semigroup |
---|
712 | of the plane curve defined by the INPUT is generated by the vectors |
---|
713 | @code{v_1,...,v_s,w_1+k*e_1,...,w_r+k*e_r}, where e_i denotes the i-th standard |
---|
714 | basis vector and k runs through all non-negative integers. The thrid entry is the conductor |
---|
715 | of the plane curve singularity. Note that r is the number of branches of the plane curve |
---|
716 | singularity and integer vectors thus have size r. |
---|
717 | NOTE: If the output is zero this means that the curve has one branch and is regular. |
---|
718 | In the reducible case the set of generators may not be minimal. |
---|
719 | @* |
---|
720 | If you are not sure whether the INPUT polynomial is reduced or not, use |
---|
721 | @code{squarefree(INPUT)} as input instead. |
---|
722 | SEE ALSO: resolutiongraph, totalmultiplicities |
---|
723 | EXAMPLE: example semigroup; shows an example |
---|
724 | " |
---|
725 | { |
---|
726 | //////////////////////////////////////////////////////////////////////////////////////// |
---|
727 | // Uses the algorithm in [CDG99] A. Campillo, F. Delgado, S.M. Gusein-Zade, On the |
---|
728 | // generators of the semigroup of a plane curve singularity, |
---|
729 | // J. London Math. Soc. (2) 60 (1999), 420-430. |
---|
730 | //////////////////////////////////////////////////////////////////////////////////////// |
---|
731 | intvec conductor; // conductor of the singularity |
---|
732 | list charexp; // characteristic exponents of the singularity |
---|
733 | int i,j; |
---|
734 | ///////////////////////////////////////////////////////////////////////////////// |
---|
735 | // Decide, which kind of input we have, and define the contact matrix |
---|
736 | ///////////////////////////////////////////////////////////////////////////////// |
---|
737 | // Input: a polynomial defining a plane curve singularity. |
---|
738 | ////////////////////////////////////////////////////////////////////////////// |
---|
739 | if (typeof(INPUT)=="poly") |
---|
740 | { |
---|
741 | /* |
---|
742 | poly FFF=squarefree(INPUT); |
---|
743 | if ( deg(FFF)!=deg(INPUT) ) |
---|
744 | { |
---|
745 | dbprint(printlevel-voice+3,"// input polynomial was not reduced"); |
---|
746 | dbprint(printlevel-voice+3,"// we continue with its reduction"); |
---|
747 | } |
---|
748 | list I@N@V=invariants(FFF); |
---|
749 | */ |
---|
750 | list I@N@V=invariants(INPUT); |
---|
751 | } |
---|
752 | else |
---|
753 | { |
---|
754 | /////////////////////////////////////////////////////////////////////////////////// |
---|
755 | // Input: the vector of characteristic exponents of an irreducible plane curve |
---|
756 | /////////////////////////////////////////////////////////////////////////////////// |
---|
757 | if (typeof(INPUT)=="intvec") |
---|
758 | { |
---|
759 | // Calculate the semigroup and the conductor directly. |
---|
760 | conductor[1]=charexp2conductor(INPUT); |
---|
761 | intvec gener=charexp2generators(INPUT); |
---|
762 | list genera; |
---|
763 | for (i=1;i<=size(gener);i++) |
---|
764 | { |
---|
765 | genera[i]=gener[i]; |
---|
766 | } |
---|
767 | return(list(genera,list(),conductor)); |
---|
768 | } |
---|
769 | else |
---|
770 | { |
---|
771 | if (typeof(INPUT)=="list") |
---|
772 | { |
---|
773 | ///////////////////////////////////////////////////////////////////////////////// |
---|
774 | // Input: intersection-matrix and characteristic exponents. |
---|
775 | ////////////////////////////////////////////////////////////////////////////// |
---|
776 | if (typeof(INPUT[1])=="intmat") |
---|
777 | { |
---|
778 | intmat contact=INPUT[1]; |
---|
779 | charexp=INPUT[2]; |
---|
780 | int n=ncols(contact); // to know how many branches we have |
---|
781 | if (n==1) // Only one branch! |
---|
782 | { |
---|
783 | return(semigroup(charexp[1])); |
---|
784 | } |
---|
785 | intmat intersecmult=charexp2inter(contact,charexp); |
---|
786 | for(i=1;i<=ncols(contact);i++) |
---|
787 | { |
---|
788 | conductor[i]=charexp2conductor(charexp[i]);//list with the characteristic exponents |
---|
789 | } |
---|
790 | } |
---|
791 | else |
---|
792 | { |
---|
793 | ///////////////////////////////////////////////////////////////////////////////// |
---|
794 | // Input: output of hnexpansion or hne in the ring created by hnexpansion |
---|
795 | ////////////////////////////////////////////////////////////////////////////// |
---|
796 | if ((typeof(INPUT[1])=="ring") or (typeof(INPUT[1])=="list")) |
---|
797 | { |
---|
798 | list I@N@V=invariants(INPUT); |
---|
799 | } |
---|
800 | else |
---|
801 | { |
---|
802 | //////////////////////////////////////////////////////////////////////////////////// |
---|
803 | // Input: output of develop or extdevelop -- irreducible plane curve singularity |
---|
804 | //////////////////////////////////////////////////////////////////////////////////// |
---|
805 | if (typeof(INPUT[1])=="matrix") |
---|
806 | { |
---|
807 | return(semigroup(invariants(INPUT)[1])); |
---|
808 | } |
---|
809 | else |
---|
810 | { |
---|
811 | ERROR("The input is invalid!"); |
---|
812 | } |
---|
813 | } |
---|
814 | } |
---|
815 | } |
---|
816 | else |
---|
817 | { |
---|
818 | ERROR("The input is invalid!"); |
---|
819 | } |
---|
820 | } |
---|
821 | } |
---|
822 | /////////////////////////////////////////////////////////////////////////////////////////////////// |
---|
823 | // If the input was a poly or an HN-Expansion, then calculate the contact matrix and char.exponents. |
---|
824 | /////////////////////////////////////////////////////////////////////////////////////////////////// |
---|
825 | if (defined(I@N@V)) |
---|
826 | { |
---|
827 | int n =size(I@N@V)-1;// number of branches |
---|
828 | // If the singularity is irreducible, then we calculate the semigroup directly. |
---|
829 | if (n==1) |
---|
830 | { |
---|
831 | return(semigroup(I@N@V[1][1])); |
---|
832 | } |
---|
833 | // If the singularity is not irreducible, then we go on. |
---|
834 | intmat contact=I@N@V[size(I@N@V)][1]; // contact matrix |
---|
835 | intmat intersecmult=I@N@V[size(I@N@V)][2]; // intersection multiplicities |
---|
836 | for(i=1;i<=n;i++) |
---|
837 | { |
---|
838 | conductor[i]=I@N@V[i][5]; |
---|
839 | charexp[i]=I@N@V[i][1]; |
---|
840 | } |
---|
841 | } |
---|
842 | ///////////////////////////////////////////////////////////////////////////////////// |
---|
843 | // If we have come so far, the curve is reducible! |
---|
844 | ///////////////////////////////////////////////////////////////////////////////////// |
---|
845 | // Compute the conductor of the curve. |
---|
846 | ///////////////////////////////////////////////////////////////////////////////////// |
---|
847 | for(i=1;i<=size(conductor);i++) |
---|
848 | { |
---|
849 | for(j=1;j<=size(conductor);j++) |
---|
850 | { |
---|
851 | conductor[i]=conductor[i]+intersecmult[i,j]; |
---|
852 | } |
---|
853 | } |
---|
854 | ///////////////////////////////////////////////////////////////////////////////////// |
---|
855 | /// The total multiplicity vectors of the exceptional divisors generate the semigroup, |
---|
856 | /// however, this system of generators is not minimal. Theorem 2 in [CDG99] leads |
---|
857 | /// to the following algorithm for minimizing the set of generators. |
---|
858 | ///////////////////////////////////////////////////////////////////////////////////// |
---|
859 | list resgr_totmult=totalmultiplicities(list(contact,charexp)); // resolution graph and tot. mult. |
---|
860 | intmat resgr=resgr_totmult[1]; // resolution graph |
---|
861 | intmat totmult=resgr_totmult[2]; // total multiplicities |
---|
862 | list conpts; // ith entry = points to which i is connected and their weights |
---|
863 | intvec series; // ith entry = row in totmult which corresponds to w_i |
---|
864 | intvec deadarc; // Star point and the point connected to the star point of a dead arc. |
---|
865 | list deadarcs; // Saves deadarc for all dead arcs. |
---|
866 | int stop,ctp,ctpstop; |
---|
867 | list ordinary_generators, series_generators; // the v_i and the w_i from the description |
---|
868 | // Find for each branch of the curve all the points in the graph to which it is connected |
---|
869 | // and save the numbers of the corresponding rows in the graph as well as the weight of the point. |
---|
870 | for (i=1;i<=ncols(resgr);i++) |
---|
871 | { |
---|
872 | conpts[i]=find_connection_points(resgr,i); |
---|
873 | } |
---|
874 | ////////////////////////////////////////////////////////////////////////////////// |
---|
875 | // Check for each point in the graph, whether it contributes to the semigroup. |
---|
876 | // If it does not contribute to the semigroup, then set the weight of the point |
---|
877 | // to zero, i.e. resgr[i,i]=0. Note that the vertex 1 always contributes! |
---|
878 | ////////////////////////////////////////////////////////////////////////////////// |
---|
879 | // Find first all the points corresponding to the branches and the end points of dead arcs. |
---|
880 | // The points to which the branch k is connected contributes to w_k. The end points of |
---|
881 | // dead arcs contribute to the v_i, while the remaining points of the dead arcs are to be removed. |
---|
882 | j=1; // Counter for dead arcs - the star points of dead arcs have to be saved for future use. |
---|
883 | for (i=2;i<=ncols(resgr);i++) |
---|
884 | { |
---|
885 | // The end points of dead arcs and the end points corresponding to the branches are |
---|
886 | // recognized by the fact that they are only connected to one other point. |
---|
887 | if (size(conpts[i][1])==1) // row i in the graph corresponds to an end point |
---|
888 | { |
---|
889 | // For an end point corresponding to a branch k the last point E_alpha(k), to which it is |
---|
890 | // connected, contributes to the generator w_k. |
---|
891 | if (resgr[i,i]<0) |
---|
892 | { |
---|
893 | series[-resgr[i,i]]=conpts[i][1][1]; // Find E_alpha(k), where k=-resgr[i,i] |
---|
894 | ctp=conpts[i][1][1]; |
---|
895 | resgr[ctp,ctp]=0; // So E_alph(k) does not contribute to the v_i. |
---|
896 | } |
---|
897 | // For an end point of a dead arc the end point contributes, but all the other points |
---|
898 | // of the dead arc - including the star point = separation pt of the dead arc - do not contribute. |
---|
899 | if (resgr[i,i]>0) |
---|
900 | { |
---|
901 | stop=0; // Stop checks whether the star point of the dead arc has been reached. |
---|
902 | ctp=conpts[i][1][1]; // The point to which the end point is connected. |
---|
903 | // Set the weights of all the other points in the dead arc to zero. |
---|
904 | while ((stop==0) and (ctp!=1)) // If the star point or the vertex 1 has been reached, stop. |
---|
905 | { |
---|
906 | deadarc[2]=i; // i might be the point connected to the star point of the dead arc. |
---|
907 | resgr[ctp,ctp]=0; |
---|
908 | if (size(conpts[ctp][1])==2) // If ctp was an ordinary vertex, we go on. |
---|
909 | { |
---|
910 | deadarc[2]=ctp; // ctp might be the point connectd to the star point of the dead arc. |
---|
911 | ctp=conpts[ctp][1][2]; // This is the second point (of higher weight) to which ctp is connected. |
---|
912 | } |
---|
913 | else // If ctp is the star point, we stop. |
---|
914 | { |
---|
915 | deadarc[1]=ctp; // Star point of this dead arc. |
---|
916 | deadarcs[j]=deadarc; // Save the j-th dead arc. |
---|
917 | j++; |
---|
918 | stop=1; |
---|
919 | } |
---|
920 | } |
---|
921 | } |
---|
922 | } |
---|
923 | } |
---|
924 | ////////////////////////////////////////////////////////////////////////////////// |
---|
925 | // Points (!=1) which are on the geodesics of every branch don't contribute. |
---|
926 | // (The geodesic of a branch is the shortest connection of the strict transform to 1 in the graph.) |
---|
927 | stop=0; // Checks whether the points on all geodesics have been found. |
---|
928 | ctp=1; // Point (=row in resgr) to be considered. |
---|
929 | int prev_ctp; // Previous point in the graph, to which ctp was connected. |
---|
930 | int dead_arc_ctp; // If ctp is the star pt of a dead arc, this is the connection pt of ctp in the d.a. |
---|
931 | int dastarptcheck; // Checks whether a point is a star point of a dead arc. |
---|
932 | if (size(conpts[1][1])>=2) // The graphs separate already at point 1. |
---|
933 | { |
---|
934 | stop=1; |
---|
935 | } |
---|
936 | else // The graphs do not separate at point 1. |
---|
937 | { |
---|
938 | prev_ctp=1; |
---|
939 | ctp=conpts[1][1][1]; // Next point to be considered. |
---|
940 | } |
---|
941 | // Pass on along the graph until we reach the first point where some branches separate. |
---|
942 | while (stop==0) |
---|
943 | { |
---|
944 | if (size(conpts[ctp][1])==2) // Point ctp is connected to 2 other points, hence is a normal vertex. |
---|
945 | { |
---|
946 | resgr[ctp,ctp]=0; // Point ctp is a normal vertex. |
---|
947 | prev_ctp=ctp; // Save the position of ctp for future use. |
---|
948 | ctp=conpts[ctp][1][2]; // Next point to which ctp is connected. |
---|
949 | } |
---|
950 | if (size(conpts[ctp][1])>3) // More than three points are connected to ctp. |
---|
951 | { |
---|
952 | resgr[ctp,ctp]=0; |
---|
953 | stop=1; // The graphs separate at point ctp. |
---|
954 | } |
---|
955 | if (size(conpts[ctp][1])==3) // At ctp a dead arc might depart or some branch(es)! |
---|
956 | { // If a dead arc departs, then the branches stay together. |
---|
957 | resgr[ctp,ctp]=0; |
---|
958 | // Check if a dead arc departs at point ctp (i.e. if ctp is the star point of a dead arc), |
---|
959 | // then the branches do not separate at ctp. |
---|
960 | dastarptcheck=0; |
---|
961 | i=1; |
---|
962 | while ((i<=size(deadarcs)) and (dastarptcheck==0)) |
---|
963 | { |
---|
964 | if (ctp==deadarcs[i][1]) // ctp is the star point of a dead arc. |
---|
965 | { |
---|
966 | dastarptcheck=1; |
---|
967 | dead_arc_ctp=deadarcs[i][2]; // The point in the dead arc to which ctp is connected. |
---|
968 | } |
---|
969 | i++; |
---|
970 | } |
---|
971 | if (dastarptcheck==0) // ctp is not the star point of a dead arc, hence the graphs separate at ctp. |
---|
972 | { |
---|
973 | stop=1; |
---|
974 | } |
---|
975 | else |
---|
976 | { |
---|
977 | // Set ctp to the point connected to ctp which is not in the dead arc and is not prev_ctp. |
---|
978 | i=1; |
---|
979 | ctpstop=0; |
---|
980 | while ((i<=3) and (ctpstop==0)) |
---|
981 | { |
---|
982 | if ((conpts[ctp][1][i]!=prev_ctp) and (conpts[ctp][1][i]!=dead_arc_ctp)) |
---|
983 | { |
---|
984 | prev_ctp=ctp; |
---|
985 | ctp=conpts[ctp][1][i]; |
---|
986 | ctpstop=1; |
---|
987 | } |
---|
988 | i++; |
---|
989 | } |
---|
990 | } |
---|
991 | } |
---|
992 | } |
---|
993 | ///////////////////////////////////////////////////////////////////////////////////////////// |
---|
994 | // Collect the generators v_j by checking which points in the graph still have |
---|
995 | // a positive weight! These points contribute their total multiplicity vector as |
---|
996 | // generator v_j. |
---|
997 | j=1; |
---|
998 | for (i=1;i<=ncols(resgr);i++) |
---|
999 | { |
---|
1000 | if (resgr[i,i]>0) |
---|
1001 | { |
---|
1002 | ordinary_generators[j]=intvec(totmult[i,1..ncols(totmult)]); |
---|
1003 | j++; |
---|
1004 | } |
---|
1005 | } |
---|
1006 | // The "exceptional" generators w_i, for which we have to include w_i+ke_i (for all k) |
---|
1007 | // are the total multiplicity vectors of the points saved in series. |
---|
1008 | for (i=1;i<=ncols(totmult);i++) |
---|
1009 | { |
---|
1010 | series_generators[i]=intvec(totmult[series[i],1..ncols(totmult)]); |
---|
1011 | } |
---|
1012 | return(list(ordinary_generators,series_generators,conductor)); |
---|
1013 | } |
---|
1014 | example |
---|
1015 | { |
---|
1016 | "EXAMPLE:"; |
---|
1017 | echo=2; |
---|
1018 | ring r=0,(x,y),ls; |
---|
1019 | // Irreducible Case |
---|
1020 | semigroup((x2-y3)^2-4x5y-x7); |
---|
1021 | // In the irreducible case, invariants() also calculates a minimal set of |
---|
1022 | // generators of the semigroup. |
---|
1023 | invariants((x2-y3)^2-4x5y-x7)[1][2]; |
---|
1024 | // Reducible Case |
---|
1025 | poly f=(y2-x3)*(y2+x3)*(y4-2x3y2-4x5y+x6-x7); |
---|
1026 | semigroup(f); |
---|
1027 | } |
---|
1028 | |
---|
1029 | |
---|
1030 | |
---|
1031 | proc charexp2generators (intvec charexp) |
---|
1032 | "USAGE: charexp2generators(v), v intvec |
---|
1033 | ASSUME: v contains the characteristic exponents of an irreducible plane |
---|
1034 | curve singularity |
---|
1035 | RETURN: intvec, the minimal set of generators of the semigroup of the plane curve singularity |
---|
1036 | SEE ALSO: invariants, resolutiongraph, totalmultiplicities, alexanderpolynomial |
---|
1037 | KEYWORDS: generators; semigroup; characteristic exponents; topological invariants |
---|
1038 | EXAMPLE: example charexp2generators; shows an example" |
---|
1039 | { |
---|
1040 | int end=size(charexp); |
---|
1041 | // If the singularity is smooth! |
---|
1042 | if (end==1) |
---|
1043 | { |
---|
1044 | return(1); |
---|
1045 | } |
---|
1046 | int i,j; |
---|
1047 | intvec gener; |
---|
1048 | intvec GGT; |
---|
1049 | for (i=1;i<=end;i++) |
---|
1050 | { |
---|
1051 | // Calculate the sequence of gcd's of the characteristic exponents. |
---|
1052 | if (i==1) |
---|
1053 | { |
---|
1054 | GGT[1]=charexp[1]; |
---|
1055 | } |
---|
1056 | else |
---|
1057 | { |
---|
1058 | GGT[i]=gcd(GGT[i-1],charexp[i]); |
---|
1059 | } |
---|
1060 | // Calculate the generators. |
---|
1061 | gener[i]=charexp[i]; |
---|
1062 | for (j=2;j<=i-1;j++) |
---|
1063 | { |
---|
1064 | gener[i]=gener[i]+((GGT[j-1]-GGT[j]) div GGT[i-1])*charexp[j]; |
---|
1065 | } |
---|
1066 | } |
---|
1067 | return(gener); |
---|
1068 | } |
---|
1069 | example |
---|
1070 | { |
---|
1071 | "EXAMPLE:"; |
---|
1072 | echo=2; |
---|
1073 | charexp2generators(intvec(28,64,66,77)); |
---|
1074 | } |
---|
1075 | |
---|
1076 | |
---|
1077 | proc charexp2multseq (intvec charexp) |
---|
1078 | "USAGE: charexp2multseq(v), v intvec |
---|
1079 | ASSUME: v contains the characteristic exponents of an irreducible plane |
---|
1080 | curve singularity |
---|
1081 | RETURN: intvec, the multiplicity sequence of the plane curve singularity |
---|
1082 | NOTE: If the curve singularity is smooth, then the multiplicity sequence is empty. |
---|
1083 | This is expressed by returning zero. |
---|
1084 | SEE ALSO: invariants, resolutiongraph, totalmultiplicities, alexanderpolynomial |
---|
1085 | KEYWORDS: characteristic exponents; multiplicity sequence; topological invariants |
---|
1086 | EXAMPLE: example charexp2multseq; shows an example" |
---|
1087 | { |
---|
1088 | int end=size(charexp); |
---|
1089 | // If the singularity is smooth! |
---|
1090 | if (end==1) |
---|
1091 | { |
---|
1092 | return(1); // ERROR: Should be 0, but for the time being, Singular returns 1. |
---|
1093 | } |
---|
1094 | intvec multseq=euclidseq(charexp[2],charexp[1]); |
---|
1095 | for (int i=3;i<=end;i++) |
---|
1096 | { |
---|
1097 | multseq=multseq,euclidseq(charexp[i]-charexp[i-1],multseq[size(multseq)]); |
---|
1098 | } |
---|
1099 | return(multseq); |
---|
1100 | } |
---|
1101 | example |
---|
1102 | { |
---|
1103 | "EXAMPLE:"; |
---|
1104 | echo=2; |
---|
1105 | charexp2multseq(intvec(28,64,66,77)); |
---|
1106 | } |
---|
1107 | |
---|
1108 | proc multseq2charexp(def v) // Procedure written by Fernando. |
---|
1109 | "USAGE: multseq2charexp(v), v intvec |
---|
1110 | ASSUME: The input is an intvec, which contains the mutiplicity sequence |
---|
1111 | of an irreducible plane curve singularity . |
---|
1112 | RETURN: An intvec, which contains the sequence of characteristic |
---|
1113 | exponents of the irreducible plane curve singularity defined by v. |
---|
1114 | EXAMPLE: example multseq2charexp; shows an example. |
---|
1115 | " |
---|
1116 | { |
---|
1117 | /////// Preamble which reduces the input of an intvec to ///////////// |
---|
1118 | /////// the originally assumed input of a list of intvecs ///////////// |
---|
1119 | /////// and tests the input. ///////////// |
---|
1120 | if (typeof(v)=="intvec") |
---|
1121 | { |
---|
1122 | list RESULT=multseq2charexp(list(v)); |
---|
1123 | return(RESULT[1]); |
---|
1124 | } |
---|
1125 | if (typeof(v)!="list") |
---|
1126 | { |
---|
1127 | ERROR("Invalid Input!"); |
---|
1128 | } |
---|
1129 | if (typeof(v)=="list") |
---|
1130 | { |
---|
1131 | int TESTV; |
---|
1132 | for (int III=1;III<=size(v);III++) |
---|
1133 | { |
---|
1134 | if (typeof(v[III])!="intvec") |
---|
1135 | { |
---|
1136 | TESTV=1; |
---|
1137 | } |
---|
1138 | } |
---|
1139 | if (TESTV==1) |
---|
1140 | { |
---|
1141 | ERROR("Invalid Input!"); |
---|
1142 | } |
---|
1143 | } |
---|
1144 | /////////////////////////////////////////////////////////// |
---|
1145 | list L=v; |
---|
1146 | int n =size(L); |
---|
1147 | // /////////////////////////////////////////////////////// |
---|
1148 | // we look the size of each vector |
---|
1149 | intvec mm; |
---|
1150 | for(int j=1;j<=n;j++) |
---|
1151 | { |
---|
1152 | mm[j]=size(L[j]); |
---|
1153 | } |
---|
1154 | // /////////////////////////////////////////////////////// |
---|
1155 | // we define some variables |
---|
1156 | list LL; |
---|
1157 | int temp, temp1,temp2; |
---|
1158 | int ind,r,l,boolean; |
---|
1159 | int old,new; |
---|
1160 | int contador; |
---|
1161 | list EUCLI,EUCLI1; |
---|
1162 | intvec exponent,exponentes1; |
---|
1163 | int new1,old1; |
---|
1164 | int contador1; |
---|
1165 | int a,b,f; |
---|
1166 | //with the for we round each branch. |
---|
1167 | for(int k=1;k<=n;k++) |
---|
1168 | { |
---|
1169 | l=1; |
---|
1170 | old=L[k][l]; |
---|
1171 | //if the vertor has more than one element |
---|
1172 | if(mm[k]<>1) |
---|
1173 | { |
---|
1174 | // /////////////////////////////////////////////////////////////////////////////// |
---|
1175 | // the first step is special because is easy to know the two first characteristic exponents |
---|
1176 | new=L[k][l+1]; |
---|
1177 | contador=1; |
---|
1178 | while(old==new)//we check how many consecutives elements are equal, starting in the first. |
---|
1179 | { |
---|
1180 | contador=contador+1; |
---|
1181 | old=new; |
---|
1182 | new=L[k][contador+1]; |
---|
1183 | } |
---|
1184 | exponent=L[k][l],contador*L[k][l]+L[k][l+contador];// those are the first two characteristic exponents. |
---|
1185 | LL[k]=exponent;// we insert them in the final matrix |
---|
1186 | EUCLI=euclides(LL[k][2],LL[k][1]);// compute the euclides algorithm for the two first characteristic exponents. |
---|
1187 | temp=size(EUCLI[1]); |
---|
1188 | // measure how many elements of the multiplicity sequence belong to the first euclidean algorithm. |
---|
1189 | for(ind=1;ind<=temp;ind=ind+1) |
---|
1190 | { |
---|
1191 | l=l+EUCLI[1][ind]; |
---|
1192 | } |
---|
1193 | l=l-1;//this is the number we are looking for. |
---|
1194 | /////////////////////////////////////////////////////////////// |
---|
1195 | r=1; |
---|
1196 | //repeat the same process until the end of the multiplicity sequence. |
---|
1197 | if(mm[k]-1>l) |
---|
1198 | { |
---|
1199 | while( l<mm[k]-1) |
---|
1200 | { |
---|
1201 | r=r+1; |
---|
1202 | old1=L[k][l]; |
---|
1203 | new1=L[k][l+1]; |
---|
1204 | contador1=0; |
---|
1205 | boolean=1; |
---|
1206 | if(old1==new1) |
---|
1207 | { |
---|
1208 | while(old1==new1 and boolean==1) |
---|
1209 | { |
---|
1210 | contador1=contador1+1; |
---|
1211 | old1=new1; |
---|
1212 | new1=L[k][l+contador1+1]; |
---|
1213 | if(size(L[k])<=l+contador1+1) |
---|
1214 | { |
---|
1215 | boolean =0; |
---|
1216 | } |
---|
1217 | } |
---|
1218 | } |
---|
1219 | temp1=size(LL[k]); |
---|
1220 | exponentes1=LL[k],LL[k][temp1]+(contador1*L[k][l])+L[k][contador1+l+1]; |
---|
1221 | LL[k]=exponentes1; |
---|
1222 | EUCLI1=euclides(LL[k][temp1+1]-LL[k][temp1],L[k][l]); |
---|
1223 | temp2=size(EUCLI1[1]); |
---|
1224 | for(ind=1;ind<=temp2;ind=ind+1) |
---|
1225 | { |
---|
1226 | l=l+EUCLI1[1][ind]; |
---|
1227 | } |
---|
1228 | } |
---|
1229 | } |
---|
1230 | } |
---|
1231 | // if the vector has only one element then the charexp is only 1. |
---|
1232 | else |
---|
1233 | { |
---|
1234 | LL[k]=1; |
---|
1235 | } |
---|
1236 | } |
---|
1237 | return(LL); |
---|
1238 | } |
---|
1239 | example |
---|
1240 | { |
---|
1241 | "EXAMPLE:";echo=2; |
---|
1242 | intvec v=2,1,1; |
---|
1243 | multseq2charexp(v); |
---|
1244 | intvec v1=4,2,2,1,1; |
---|
1245 | multseq2charexp(v1); |
---|
1246 | } |
---|
1247 | |
---|
1248 | proc charexp2inter (intmat contact, list charexp) |
---|
1249 | "USAGE: charexp2inter(contact,charexp), contact matrix, charexp list |
---|
1250 | ASSUME: charexp contains the integer vectors of characteristic exponents |
---|
1251 | of the branches of a plane curve singularity, and contact is their |
---|
1252 | contact matrix |
---|
1253 | RETURN: intmat, the matrix intersection multiplicities of the branches |
---|
1254 | SEE ALSO: invariants, resolutiongraph, totalmultiplicities, semigroup |
---|
1255 | KEYWORDS: contact matrix; characteristic exponents; intersection multiplicity; topological invariants |
---|
1256 | EXAMPLE: example charexp2inter; shows an example" |
---|
1257 | { |
---|
1258 | int n=ncols(contact); |
---|
1259 | int i,j,k; |
---|
1260 | list multpl; |
---|
1261 | int max=0; |
---|
1262 | intvec helpvect; |
---|
1263 | intmat inters[n][n]; |
---|
1264 | // Calculate the multiplicity sequences of the branches. |
---|
1265 | for (i=1;i<=n;i++) |
---|
1266 | { |
---|
1267 | multpl[i]=charexp2multseq(charexp[i]); |
---|
1268 | // Find the maximal length of a multiplicity sequence. |
---|
1269 | if (max<size(multpl[i])) |
---|
1270 | { |
---|
1271 | max=size(multpl[i]); |
---|
1272 | } |
---|
1273 | } |
---|
1274 | // If the contact of certain branches is higher than the maximal length of the |
---|
1275 | // multiplicity sequence, then max should be the maximal contact! |
---|
1276 | helpvect=max,intvec(contact); |
---|
1277 | max=max_in_intvec(helpvect)[1]; |
---|
1278 | // Prolong them by 1s, in order to take care of higher contact. |
---|
1279 | for (i=1;i<=n;i++) |
---|
1280 | { |
---|
1281 | helpvect=multpl[i]; |
---|
1282 | for (j=size(multpl[i]+1);j<=max;j++) |
---|
1283 | { |
---|
1284 | helpvect[j]=1; |
---|
1285 | } |
---|
1286 | multpl[i]=helpvect; |
---|
1287 | } |
---|
1288 | // Calculate the intersection numbers of the branches: for two branches f_i and f_j |
---|
1289 | // this is the sum over mult_k(f_i)*mult_k(f_j), where k runs over all infinitely |
---|
1290 | // near points which f_i and f_j share, i.e. from 1 to contact[i,j]. |
---|
1291 | for (i=1;i<=n;i++) |
---|
1292 | { |
---|
1293 | for (j=i+1;j<=n;j++) |
---|
1294 | { |
---|
1295 | for (k=1;k<=contact[i,j];k++) |
---|
1296 | { |
---|
1297 | inters[i,j]=inters[i,j]+multpl[i][k]*multpl[j][k]; |
---|
1298 | } |
---|
1299 | inters[j,i]=inters[i,j]; |
---|
1300 | } |
---|
1301 | } |
---|
1302 | return(inters); |
---|
1303 | } |
---|
1304 | example |
---|
1305 | { |
---|
1306 | "EXAMPLE:";echo=2; |
---|
1307 | ring r=0,(x,y),ds; |
---|
1308 | list INV=invariants((x2-y3)*(x3-y2)*((x2-y3)^2-4x5y-x7)); |
---|
1309 | intmat contact=INV[4][1]; |
---|
1310 | list charexp=INV[1][1],INV[2][1],INV[3][1]; |
---|
1311 | // The intersection matrix is INV[4][2]. |
---|
1312 | print(INV[4][2]); |
---|
1313 | // And it is calulated as ... |
---|
1314 | print(charexp2inter(contact,charexp)); |
---|
1315 | } |
---|
1316 | |
---|
1317 | |
---|
1318 | proc charexp2conductor(intvec B) // Procedure written by Fernando |
---|
1319 | "USAGE: charexp2conductor(v), v intvec |
---|
1320 | ASSUME: v contains the characteristic exponents of an irreducible plane |
---|
1321 | curve singularity |
---|
1322 | RETURN: int, the conductor of the plane curve singularity |
---|
1323 | NOTE: If the curve singularity is smooth, the conductor is zero. |
---|
1324 | SEE ALSO: invariants, resolutiongraph, totalmultiplicities, semigroup |
---|
1325 | KEYWORDS: conductor; characteristic exponents; multiplicity sequence; topological invariants |
---|
1326 | EXAMPLE: example charexp2conductor; shows an example" |
---|
1327 | { |
---|
1328 | intvec E; |
---|
1329 | int i,conductor; |
---|
1330 | E[1]=B[1]; |
---|
1331 | for(i=2;i<=size(B);i++) |
---|
1332 | { |
---|
1333 | E[i]=gcd(E[i-1],B[i]); |
---|
1334 | } |
---|
1335 | conductor=1-E[1]; |
---|
1336 | for(i=2;i<=size(B);i++) |
---|
1337 | { |
---|
1338 | conductor=conductor+(B[i]*(E[i-1]-E[i])); |
---|
1339 | } |
---|
1340 | return(conductor); |
---|
1341 | } |
---|
1342 | example |
---|
1343 | { |
---|
1344 | "EXAMPLE:"; |
---|
1345 | echo=2; |
---|
1346 | charexp2conductor(intvec(2,3)); // A1-Singularity |
---|
1347 | charexp2conductor(intvec(28,64,66,77)); |
---|
1348 | } |
---|
1349 | |
---|
1350 | |
---|
1351 | proc charexp2poly(intvec v, vector a) // Procedure written by Fernando. |
---|
1352 | "USAGE: charexp2poly(v,a); intvec v, vector a. |
---|
1353 | ASSUME: v an intvec containing the characterictic exponents of an irreducible plane curve singularity. |
---|
1354 | a a vector containing the coefficients of a parametrization given by x(t)=x^v[1], |
---|
1355 | y(t)=a(1)t^v[2]+...+a[n-1]t^v[n], i.e. the entries of a are of type number. |
---|
1356 | RETURN: A polynomial f in the first two variables of the basering, such that f defines an |
---|
1357 | irreducible plane curve singularity with characteristic exponents v. |
---|
1358 | NOTE: The entries in a should be of type number and the vector v should be the sequence of |
---|
1359 | characteristic exponents of an irreducible plane curve singularity in order to |
---|
1360 | get a sensible result, |
---|
1361 | SEE ALSO: charexp2multseq, multseq2charexp. |
---|
1362 | EXAMPLE: example charexp2poly; shows an example |
---|
1363 | " |
---|
1364 | { |
---|
1365 | int n=size(v); |
---|
1366 | vector expo; |
---|
1367 | int i,j,s; |
---|
1368 | for(i=1;i<=v[1];i++) |
---|
1369 | { |
---|
1370 | expo[i]=0;//initialize to 0. |
---|
1371 | } |
---|
1372 | for(i=2;i<=n;i++) |
---|
1373 | { |
---|
1374 | s=v[i] mod v[1];//calculate the position. |
---|
1375 | expo=expo-a[i-1]*var(1)^((v[i]-s) div v[1])*gen(s+1);//save in expo -var(1) to the power the corresponding |
---|
1376 | //but only in the right positions |
---|
1377 | } |
---|
1378 | matrix M[v[1]][v[1]]; |
---|
1379 | //construct the matrix that generates the polynomial |
---|
1380 | for(i=1;i<=v[1];i++) |
---|
1381 | { |
---|
1382 | M[i,i]=var(2)+expo[1];//The diagonal |
---|
1383 | for(j=1;j<=v[1]-i;j++) |
---|
1384 | { |
---|
1385 | M[j,j+i]=expo[i+1];//over diagonal |
---|
1386 | } |
---|
1387 | for(j=1;j<=v[1]-i;j++) |
---|
1388 | { |
---|
1389 | M[j+i,j]=var(1)*expo[1+v[1]-i];//under diagonal |
---|
1390 | } |
---|
1391 | } |
---|
1392 | //the poynomial is the determinant of the matrix |
---|
1393 | poly irredpoly=det(M); |
---|
1394 | return(irredpoly) |
---|
1395 | } |
---|
1396 | example |
---|
1397 | { |
---|
1398 | "EXAMPLE:";echo=2; |
---|
1399 | ring r=0,(x,y),dp; |
---|
1400 | intvec v=8,12,14,17; |
---|
1401 | vector a=[1,1,1]; |
---|
1402 | poly f=charexp2poly(v,a); |
---|
1403 | f; |
---|
1404 | invariants(f)[1][1]; // The characteristic exponents of f. |
---|
1405 | } |
---|
1406 | |
---|
1407 | proc tau_es2 (def INPUT, list #) |
---|
1408 | "USAGE: tau_es2(INPUT); INPUT poly or list |
---|
1409 | ASSUME: INPUT is either a REDUCED bivariate polynomial defining a plane curve singularity, |
---|
1410 | or the output of @code{hnexpansion(f[,\"ess\"])}, or the list @code{hne} in |
---|
1411 | the ring created by @code{hnexpansion(f[,\"ess\"])}, or the output of |
---|
1412 | @code{develop(f)} resp. of @code{extdevelop(f,n)}, or a list containing |
---|
1413 | the contact matrix and a list of integer vectors with the characteristic exponents |
---|
1414 | of the branches of a plane curve singularity, or an integer vector containing the |
---|
1415 | characteristic exponents of an irreducible plane curve singularity. |
---|
1416 | RETURN: int, the equisingular Tjurina number of f, i. e. the codimension of the mu-constant |
---|
1417 | stratum in the semiuniversal deformation of f, where mu is the Milnor number of f. |
---|
1418 | NOTE: The equisingular Tjurina number is calculated with the aid of a Hamburger-Noether |
---|
1419 | expansion, which is the hard part of the calculation. |
---|
1420 | In case the Hamburger-Noether expansion of the curve f is needed |
---|
1421 | for other purposes as well it is better to calculate this first |
---|
1422 | with the aid of @code{hnexpansion} and use it as input instead of |
---|
1423 | the polynomial itself. |
---|
1424 | @* |
---|
1425 | If you are not sure whether the INPUT polynomial is reduced or not, use |
---|
1426 | @code{squarefree(INPUT)} as input instead. |
---|
1427 | SEE ALSO: tau_es, develop, hnexpansion, totalmultiplicities, equising_lib |
---|
1428 | EXAMPLE: example tau_es2; shows an example |
---|
1429 | " |
---|
1430 | { |
---|
1431 | // If the input is a weighted homogeneous polynomial, then use a direct algorithm to |
---|
1432 | // calculate the equisingular Tjurina number, by caluclating a K-basis of the |
---|
1433 | // Tjurina algebra and omitting those elements with weighted degree at least the |
---|
1434 | // weighted degree of the polynomial. -- If an additional input # is given (which is |
---|
1435 | // not just the number 1 !!!), the procedure always uses the recursive algorithm. |
---|
1436 | if ((typeof(INPUT)=="poly") and (size(#)==0)) |
---|
1437 | { |
---|
1438 | if (qhweight(INPUT)[1]!=0) |
---|
1439 | { |
---|
1440 | /* |
---|
1441 | poly f=squarefree(INPUT); |
---|
1442 | if ( deg(f)!=deg(INPUT) ) |
---|
1443 | { |
---|
1444 | dbprint(printlevel-voice+3,"// input polynomial was not reduced"); |
---|
1445 | dbprint(printlevel-voice+3,"// we continue with its reduction"); |
---|
1446 | } |
---|
1447 | return(tau_es_qh(f)); |
---|
1448 | */ |
---|
1449 | return(tau_es_qh(INPUT)); |
---|
1450 | } |
---|
1451 | } |
---|
1452 | // Else apply the recursive algorithm from Eugenii Shustin, On Manifolds of Singular |
---|
1453 | // Algebraic Curves, Selecta Math. Sov. Vol 10, No. 1 (1991), p. 31. |
---|
1454 | int i,j,k; // Laufvariablen |
---|
1455 | int tau,tau_i; // Variable for es-Tjurina no., multiplicitiy, and es-Tjurina at blow ups |
---|
1456 | intmat contact; // contact matrix |
---|
1457 | list graphs, multipl; // resolution graphs and multiplicity sequences |
---|
1458 | list graphs_ed, multipl_ed; // res.graphs and mult.seq. of one curve on 1st ex. div. |
---|
1459 | intmat contact_ed; // contact matrix of one curve on 1st exceptional divisor |
---|
1460 | int d_i; // correction term |
---|
1461 | intvec connections; |
---|
1462 | intmat graph; |
---|
1463 | intvec multseq; |
---|
1464 | int s,e; |
---|
1465 | int multi; // multiplicity of the curve defined by the input. |
---|
1466 | ///////////////////////////////////////////////////////// |
---|
1467 | // Check what type the input is, and act accordingly. |
---|
1468 | if (typeof(INPUT)=="list") |
---|
1469 | { |
---|
1470 | if (size(INPUT)==3) |
---|
1471 | { |
---|
1472 | // If the INPUT is the output of totalmultiplicities(INPUT,"tau") |
---|
1473 | if ((typeof(INPUT[1])=="list") and (typeof(INPUT[2])=="list") and (typeof(INPUT[3])=="intmat")) |
---|
1474 | { |
---|
1475 | graphs=INPUT[1]; |
---|
1476 | multipl=INPUT[2]; |
---|
1477 | contact=INPUT[3]; |
---|
1478 | } |
---|
1479 | // Otherwise call the procedure with the output of totalmultiplicities(INPUT,"tau"). |
---|
1480 | else |
---|
1481 | { |
---|
1482 | return(tau_es2(totalmultiplicities(INPUT,"tau"))); |
---|
1483 | } |
---|
1484 | } |
---|
1485 | else |
---|
1486 | { |
---|
1487 | return(tau_es2(totalmultiplicities(INPUT,"tau"))); |
---|
1488 | } |
---|
1489 | } |
---|
1490 | // Otherwise call the procedure with the output of totalmultiplicities(INPUT,"tau"). |
---|
1491 | else |
---|
1492 | { |
---|
1493 | return(tau_es2(totalmultiplicities(INPUT,"tau"))); |
---|
1494 | } |
---|
1495 | /// End of checking the input |
---|
1496 | /////////////////////////////////////////////////////// |
---|
1497 | // If the singularity is smooth, the equisingular Tjurina number is zero. |
---|
1498 | if ((ncols(contact)==1) and (multipl[1][1]<=1)) |
---|
1499 | { |
---|
1500 | return(0); |
---|
1501 | } |
---|
1502 | // Otherwise calculate the multiplicity of the singularity. |
---|
1503 | for (i=1;i<=size(multipl);i++) |
---|
1504 | { |
---|
1505 | multi=multi+multipl[i][1]; |
---|
1506 | } |
---|
1507 | // Find the branches which stay together after blowing up once, and group them together. |
---|
1508 | k=0; |
---|
1509 | intvec curves_on_first_ex_div=k; // If this is 0=i_0,i_1,...i_s=ncols(contact), then the branches |
---|
1510 | while (k<ncols(contact)) // (i_j)+1,...,i_(j+1) stay together after the first blowing up, |
---|
1511 | { // and s is the number of infinitely near points of first order. |
---|
1512 | k=find_last_non_one(intvec(contact[k+1,1..ncols(contact)]),k+2); |
---|
1513 | curves_on_first_ex_div=curves_on_first_ex_div,k; |
---|
1514 | } |
---|
1515 | // And then calculate the equisingular Tjurina numbers in the infinitely near points of |
---|
1516 | // first order plus some correction term (= the number of blow ups needed to separate the |
---|
1517 | // strict transform from the exceptional divisor in the corresponding infinitely near point) |
---|
1518 | // and add them. |
---|
1519 | for (i=1;i<=size(curves_on_first_ex_div)-1;i++) |
---|
1520 | { |
---|
1521 | s=curves_on_first_ex_div[i]+1; |
---|
1522 | e=curves_on_first_ex_div[i+1]; |
---|
1523 | graphs_ed=list(graphs[s..e]); |
---|
1524 | multipl_ed=list(multipl[s..e]); |
---|
1525 | contact_ed=intmat_minus_one(intmat(intvec(contact[s..e,s..e]),e-s+1,e-s+1)); |
---|
1526 | connections=0; |
---|
1527 | // Adjust the graphs and multiplicity sequences by cutting the first level. |
---|
1528 | // And find in each graph the point to which the point 1 is connected = 1 + number of |
---|
1529 | // blow ups needed to separate the strict transform from the first exceptional divisor. |
---|
1530 | for (j=1;j<=size(graphs_ed);j++) |
---|
1531 | { |
---|
1532 | // Adjust the graphs and find the connection points. |
---|
1533 | graph=graphs_ed[j]; |
---|
1534 | connections[j]=find_connection_point(intvec(graph[1,1..ncols(graph)]),1); |
---|
1535 | graph=delete_row(delete_col(graph,1),1); |
---|
1536 | graphs_ed[j]=graph; |
---|
1537 | // Adjust the multiplicity sequences. |
---|
1538 | multseq=multipl_ed[j]; |
---|
1539 | if (size(multseq)==1) // If multseq has only one entry, then their is only |
---|
1540 | { // one branch left and it is smooth! So set multipl_ed[j]=1. |
---|
1541 | multipl_ed[j]=intvec(1); |
---|
1542 | } |
---|
1543 | else // Otherwise just cut the first entry. |
---|
1544 | { |
---|
1545 | multipl_ed[j]=intvec(multseq[2..size(multseq)]); |
---|
1546 | } |
---|
1547 | } |
---|
1548 | // Calculate the equisingular Tjurina number of the strict transform at the i-th |
---|
1549 | // intersection point with the first exceptional divisor. |
---|
1550 | tau_i=tau_es2(list(graphs_ed,multipl_ed,contact_ed)); |
---|
1551 | // Calculate the number no. of blow ups needed to separate the strict transform of |
---|
1552 | // the curve defined by the input from the first exceptional divisor at their i-th |
---|
1553 | // intersection point. |
---|
1554 | d_i=max_in_intvec(connections)-1; |
---|
1555 | // If d_i is negative, then there was only one branch left and it was smooth, so d_i should be 1. |
---|
1556 | if (d_i<0) |
---|
1557 | { |
---|
1558 | d_i=1; |
---|
1559 | } |
---|
1560 | // If the curve defined by graphs_ed was smooth, then the summand has to be reduced by 1. |
---|
1561 | if (tau_i==0) |
---|
1562 | { |
---|
1563 | tau_i=-1; |
---|
1564 | } |
---|
1565 | tau=tau+tau_i+d_i; |
---|
1566 | } |
---|
1567 | // The equisingular Tjurina number is then calculated by adding the following term. |
---|
1568 | tau=tau+((multi*(multi+1)) div 2)-2; |
---|
1569 | return(tau); |
---|
1570 | } |
---|
1571 | example |
---|
1572 | { |
---|
1573 | "EXAMPLE:"; |
---|
1574 | echo=2; |
---|
1575 | ring r=0,(x,y),ls; |
---|
1576 | poly f1=y2-x3; |
---|
1577 | poly f2=(y2-x3)^2-4x5y-x7; |
---|
1578 | poly f3=y3-x2; |
---|
1579 | tau_es2(f1); |
---|
1580 | tau_es2(f2); |
---|
1581 | tau_es2(f1*f2*f3); |
---|
1582 | } |
---|
1583 | |
---|
1584 | |
---|
1585 | /////////////////////////////////////////////////////////////////////////////////////////// |
---|
1586 | // Static procedures. |
---|
1587 | /////////////////////////////////////////////////////////////////////////////////////////// |
---|
1588 | |
---|
1589 | static proc euclidseq(int a,int b) |
---|
1590 | "USAGE: euclidseq(a,b), a,b integers |
---|
1591 | RETURN: intvec, the divisors in the euclidean alogrithm with multiplicities |
---|
1592 | KEYWORDS: Euclidean algorithm; multiplicity sequence |
---|
1593 | NOTE : This procedure is for internal use only; it is called by charexp2multseq. |
---|
1594 | " |
---|
1595 | { |
---|
1596 | int i; |
---|
1597 | // multseq saves in each step of the Euclidean algorithm q-times the divisor b |
---|
1598 | intvec multseq; |
---|
1599 | int q=a div b; |
---|
1600 | int r=a mod b; |
---|
1601 | for (i=1;i<=q;i++) |
---|
1602 | { |
---|
1603 | multseq[i]=b; |
---|
1604 | } |
---|
1605 | int s=q; // size of multseq |
---|
1606 | a=b; |
---|
1607 | b=r; |
---|
1608 | while(r<>0) |
---|
1609 | { |
---|
1610 | q=a div b; |
---|
1611 | r=a mod b; |
---|
1612 | for (i=1;i<=q;i++) |
---|
1613 | { |
---|
1614 | multseq[s+i]=b; |
---|
1615 | } |
---|
1616 | s=s+q; // size of multseq |
---|
1617 | a=b; |
---|
1618 | b=r; |
---|
1619 | } |
---|
1620 | return(multseq); |
---|
1621 | } |
---|
1622 | |
---|
1623 | static proc puiseuxchainpart (int piij, int muij, intvec multpl, int initial_tm, int end_tm, int j) |
---|
1624 | "USAGE: puiseuxchainpart(piij,muiij,multpl,initial_tm,end_tm,j); |
---|
1625 | RETURN: list L, L[1] incidence matrix of a part of the Puiseux chain, L[2] the total |
---|
1626 | multiplicities of this part of the Puiseux chain. |
---|
1627 | NOTE: This procedure is only for internal use; it is called by puiseuxchain. |
---|
1628 | " |
---|
1629 | { |
---|
1630 | int delta=1; |
---|
1631 | if (j==1){delta=0;} // Delta measures whether j is 1 or not. |
---|
1632 | intvec totalmultiplicity; // Keeps the total multiplicities. |
---|
1633 | intmat pcp[muij][muij]; // Keeps the incidence matrix of the Puiseuxchainpart S_i,j. |
---|
1634 | // Calculate the total multiplicities and the Puiseuxchainpart S_i,j. |
---|
1635 | totalmultiplicity[1]=initial_tm+end_tm+multpl[1]; |
---|
1636 | pcp[1,1]=piij+1; |
---|
1637 | for (int k=2;k<=muij;k++) |
---|
1638 | { |
---|
1639 | pcp[k,k]=piij+k; |
---|
1640 | pcp[k-1,k]=1; |
---|
1641 | pcp[k,k-1]=1; |
---|
1642 | totalmultiplicity[k]=totalmultiplicity[k-1]+delta*initial_tm+multpl[k]; |
---|
1643 | } |
---|
1644 | list result=pcp,totalmultiplicity; |
---|
1645 | return(result); |
---|
1646 | } |
---|
1647 | |
---|
1648 | static proc puiseuxchain (int initial, intvec divseq, intvec multpl, int initial_tm) |
---|
1649 | "USAGE: puiseuxchain(initial,divseq,multpl,initial_tm); int initial, initial_tm, intvec divseq, multpl |
---|
1650 | RETURN: list L, L[1] incidence matrix of a Puiseux chain, L[2] the weight of the point to which the |
---|
1651 | previous Puiseux chain has to be connected, L[3] the sequence of total multiplicities of |
---|
1652 | the points in this Puiseux chain. |
---|
1653 | NOTE: This procedure is only for internal use; it is called by irred_resgraph_totmult. |
---|
1654 | " |
---|
1655 | { |
---|
1656 | int j,k,l,connectpoint; |
---|
1657 | intvec multpli; |
---|
1658 | int pc_tm=initial_tm; // Keeps the total multipl. of the endpoint of P_i-1. |
---|
1659 | int end_tm=0; |
---|
1660 | int start=1; |
---|
1661 | int omega=size(divseq); |
---|
1662 | // Keep the endpoints of the puiseuxchainparts (minus initial) s_i,j in divseq_sum. |
---|
1663 | intvec divseq_sum=divseq[1]; |
---|
1664 | for (j=2;j<=omega;j++) |
---|
1665 | { |
---|
1666 | divseq_sum[j]=divseq_sum[j-1]+divseq[j]; |
---|
1667 | } |
---|
1668 | // Define the connecting point of the Puiseuxchain P_i-1 with P_i. |
---|
1669 | if (divseq[1]==0) |
---|
1670 | { |
---|
1671 | // If divseq[1]=mu_i,1=0, then the Puiseuxchainpart S_i,1 is empty. |
---|
1672 | // We may start building the Puiseuxchain with part 2, i.e. set start=2. |
---|
1673 | start=2; |
---|
1674 | if (omega>=3) |
---|
1675 | { |
---|
1676 | connectpoint=initial+divseq_sum[2]+1; // startpoint of s_i+1,3 |
---|
1677 | } |
---|
1678 | else |
---|
1679 | { |
---|
1680 | connectpoint=initial+divseq_sum[2]; // endpoint of s_i+1,2 |
---|
1681 | } |
---|
1682 | } |
---|
1683 | else |
---|
1684 | { |
---|
1685 | connectpoint=initial+1; |
---|
1686 | } |
---|
1687 | // Build the Puiseuxchainparts s_i,j and put them as blocks into pc=P_i. |
---|
1688 | multpli=multpl[initial+1..initial+divseq_sum[start]]; |
---|
1689 | list pcp=puiseuxchainpart(initial,divseq[start],multpli,initial_tm,end_tm,start); |
---|
1690 | intmat pc=pcp[1]; |
---|
1691 | intvec tm=pcp[2]; |
---|
1692 | for (j=start+1;j<=omega;j++) |
---|
1693 | { |
---|
1694 | // Keep the final total multipl. of the puiseuxchainpart S_i,j-2 resp. P_i-1, if S_i,1 empty. |
---|
1695 | if (j>2){end_tm=initial_tm;} |
---|
1696 | // Calculate the endpoint of S_i,j-1. |
---|
1697 | initial=initial+divseq[j-1]; |
---|
1698 | // Keep the total multiplicity of the endpoint of S_i,j-1 |
---|
1699 | initial_tm=pcp[2][size(pcp[2])]; |
---|
1700 | // Build the new puiseuxchainpart S_i,j |
---|
1701 | multpli=multpl[initial+1..initial+divseq[j]]; |
---|
1702 | pcp=puiseuxchainpart(initial,divseq[j],multpli,initial_tm,end_tm,j); |
---|
1703 | pc=addmat(pc,pcp[1]); |
---|
1704 | tm=tm,pcp[2]; |
---|
1705 | } |
---|
1706 | // Connect the Puiseuxchainparts s_i,j. |
---|
1707 | for (j=start;j<=omega-2;j++) |
---|
1708 | { |
---|
1709 | k=divseq_sum[j]; // endpoint of s_i,j |
---|
1710 | l=divseq_sum[j+1]+1; // startpoint of s_i,j+2 |
---|
1711 | pc[k,l]=1; // connecting these |
---|
1712 | pc[l,k]=1; |
---|
1713 | } |
---|
1714 | if (omega>=start+1) // If there are at least two non-empty s_i,j. |
---|
1715 | { |
---|
1716 | k=divseq_sum[omega-1]; // endpoint of s_i,omega-1 |
---|
1717 | l=divseq_sum[omega]; // endpoint of s_i,omega |
---|
1718 | pc[k,l]=1; // connecting these |
---|
1719 | pc[l,k]=1; |
---|
1720 | } |
---|
1721 | list ergebnis=pc,connectpoint,tm; |
---|
1722 | return(ergebnis); |
---|
1723 | } |
---|
1724 | |
---|
1725 | static proc irred_resgraph_totmult (intvec charexp) |
---|
1726 | "USAGE: irred_resgraph_totmult(charexp); charexp intvec |
---|
1727 | ASSUME: charexp is an integer vector containg the characteristic exponents of |
---|
1728 | an irreducible plane curve singularity |
---|
1729 | RETURN: list L, L[1] is the incidence matrix of the resolution graph of the plane curve |
---|
1730 | singularity defined by INPUT, and L[2] is its sequence of total multiplicities |
---|
1731 | NOTE: This procedure is only for internal use; it is called by resgraph. |
---|
1732 | " |
---|
1733 | { |
---|
1734 | int k,l; |
---|
1735 | intvec multpl=charexp2multseq(charexp); // multiplicity sequence of the singularity |
---|
1736 | // Do first the case where the singularity is actually smooth. |
---|
1737 | if (size(charexp)==1) |
---|
1738 | { |
---|
1739 | intmat resgraph[1][1]=0; |
---|
1740 | intvec tm=1; // there is no exceptional divisor in the resolution - ERROR: should be 0, |
---|
1741 | // but for the time being, Singular returns 1 as multiplicity of smooth curves |
---|
1742 | list result=resgraph,tm,multpl; |
---|
1743 | return(result); |
---|
1744 | } |
---|
1745 | // End of the smooth case |
---|
1746 | int initial_tm=0; // total multipl. of the endpoint of Puiseux chain P_i-1 |
---|
1747 | int g=size(charexp); |
---|
1748 | list es=divsequence(charexp[2],charexp[1]); // keeps the lengths of the Puiseuxchainparts s_i,j |
---|
1749 | intvec divseq=es[1]; |
---|
1750 | int r=es[2]; |
---|
1751 | int initial=0; |
---|
1752 | // Build the Puiseuxchains P_i and put them as blocks into a matrix. |
---|
1753 | list pc=puiseuxchain(initial,divseq,multpl,initial_tm); |
---|
1754 | intmat resgraph=pc[1]; |
---|
1755 | intvec endpoints=resgraph[nrows(resgraph),ncols(resgraph)]; |
---|
1756 | intvec connectpoints=pc[2]; |
---|
1757 | intvec tm=pc[3]; |
---|
1758 | for (int i=3;i<=g;i++) |
---|
1759 | { |
---|
1760 | initial_tm=tm[size(tm)]; |
---|
1761 | es=divsequence(charexp[i]-charexp[i-1],r); |
---|
1762 | divseq=es[1]; |
---|
1763 | r=es[2]; |
---|
1764 | initial=endpoints[size(endpoints)]; |
---|
1765 | pc=puiseuxchain(initial,divseq,multpl,initial_tm); |
---|
1766 | resgraph=addmat(resgraph,pc[1]); |
---|
1767 | endpoints=endpoints,resgraph[nrows(resgraph),ncols(resgraph)]; |
---|
1768 | connectpoints=connectpoints,pc[2]; |
---|
1769 | tm=tm,pc[3]; |
---|
1770 | } |
---|
1771 | // Adding the * for the strict transform to the Graph. |
---|
1772 | resgraph=addmat(resgraph,0); |
---|
1773 | // The connecting point of the * with the graph. |
---|
1774 | connectpoints=connectpoints,nrows(resgraph); |
---|
1775 | // Connect the P_i with each other and with *. |
---|
1776 | for (i=2;i<=g;i++) |
---|
1777 | { |
---|
1778 | k=endpoints[i-1]; // endpoint of P_i-1 |
---|
1779 | l=connectpoints[i]; // conncting point of P_i resp. of * |
---|
1780 | resgraph[k,l]=1; // connecting these |
---|
1781 | resgraph[l,k]=1; |
---|
1782 | } |
---|
1783 | list result=resgraph,tm,multpl; //HIER GEAENDERT!!!! |
---|
1784 | return(result); |
---|
1785 | } |
---|
1786 | |
---|
1787 | |
---|
1788 | static proc max_in_intvec (intvec v, list #) |
---|
1789 | "USAGE: max_in_intvec(v); v intvec |
---|
1790 | RETURN: int m, maximum of the integers in v |
---|
1791 | USAGE: max_in_intvec(v,1); v intvec |
---|
1792 | RETURN: intvec m, m[1] maximum of the integers in v, m[2] position of the |
---|
1793 | last occurence of the maximum in v |
---|
1794 | NOTE: This procedure is only for internal use; this procedure is called by |
---|
1795 | totalmultiplicities and semigroup. |
---|
1796 | " |
---|
1797 | { |
---|
1798 | int max=v[1]; |
---|
1799 | int maxpos=1; |
---|
1800 | for (int i=2;i<=size(v);i++) |
---|
1801 | { |
---|
1802 | if (v[i]>max) |
---|
1803 | { |
---|
1804 | max=v[i]; |
---|
1805 | maxpos=i; |
---|
1806 | } |
---|
1807 | } |
---|
1808 | if (size(#)==0) |
---|
1809 | { |
---|
1810 | return(max); |
---|
1811 | } |
---|
1812 | else |
---|
1813 | { |
---|
1814 | return(intvec(max,maxpos)); |
---|
1815 | } |
---|
1816 | } |
---|
1817 | |
---|
1818 | static proc addmat (intmat A,intmat B) |
---|
1819 | "USAGE: max_in_intvec(A,B); A, B integer matrices |
---|
1820 | RETURN: intmat C, block matrix with left-upper block A, right-lower block B |
---|
1821 | NOTE: This procedure is only for internal use; this procedure is called several times. |
---|
1822 | " |
---|
1823 | { |
---|
1824 | int nc=ncols(A); |
---|
1825 | int nr=nrows(A); |
---|
1826 | int mc=ncols(B); |
---|
1827 | int mr=nrows(B); |
---|
1828 | int i,j; |
---|
1829 | intmat AB[nr+mr][nc+mc]; |
---|
1830 | for (i=1;i<=nr;i++) |
---|
1831 | { |
---|
1832 | for (j=1;j<=nc;j++) |
---|
1833 | { |
---|
1834 | AB[i,j]=A[i,j]; |
---|
1835 | } |
---|
1836 | } |
---|
1837 | for (i=1;i<=mr;i++) |
---|
1838 | { |
---|
1839 | for (j=1;j<=mc;j++) |
---|
1840 | { |
---|
1841 | AB[i+nr,j+nc]=B[i,j]; |
---|
1842 | } |
---|
1843 | } |
---|
1844 | return(AB); |
---|
1845 | } |
---|
1846 | |
---|
1847 | static proc divsequence(int a,int b) |
---|
1848 | "USAGE: divsequence(a,b); a,b integers |
---|
1849 | RETURN: list l, l[1] the multiplicities of the divisors in the Euclidean algorithm |
---|
1850 | and l[2] the last non-zero remainder in the Euclidean algorithm |
---|
1851 | NOTE: This procedure is only for internal use; it is called in irred_res_graph. |
---|
1852 | " |
---|
1853 | { |
---|
1854 | int q=a div b; |
---|
1855 | int r=a mod b; |
---|
1856 | intvec divseq=q; |
---|
1857 | while(r<>0) |
---|
1858 | { |
---|
1859 | a=b; |
---|
1860 | b=r; |
---|
1861 | q=a div b; |
---|
1862 | r=a mod b; |
---|
1863 | divseq = divseq,q; |
---|
1864 | } |
---|
1865 | list result=divseq,b; |
---|
1866 | return(result); |
---|
1867 | } |
---|
1868 | |
---|
1869 | |
---|
1870 | |
---|
1871 | static proc adjust_tot_mult (intvec rtm_fix, rtm_var, rmt_fix, rmt_var, p, q, stricttransforms, int k) |
---|
1872 | "USAGE: adjust_tot_mult(v1,v2,v3,v4,p,q,st,k); v1,...,st intvecs and k an integer |
---|
1873 | RETURN: intvec rtm_var, adjusted total multiplicities |
---|
1874 | NOTE: This procedure is only for internal use; it is called in totalmultiplicities. |
---|
1875 | " |
---|
1876 | { |
---|
1877 | int j,l,store; // Help variables. |
---|
1878 | // Recalculate the entries in rtm_var from stricttransforms[k]+1,...,stricttransforms[k+1]-1. |
---|
1879 | for (j=stricttransforms[k]+1;j<stricttransforms[k+1];j++) |
---|
1880 | { |
---|
1881 | if (rtm_var[j]==0) // If the entry is non-zero, we know that it is already correct. |
---|
1882 | { |
---|
1883 | // Check if the vertex in the fixed part is connected to one or to two vertices of lower weight. |
---|
1884 | if (j==stricttransforms[k]+1) // The vertex of weight 1 less is p[k], to which the subgraph is glued. |
---|
1885 | { |
---|
1886 | store=rtm_fix[j]-rmt_fix[j]-rtm_fix[p[k]]; |
---|
1887 | } |
---|
1888 | else // The vertex of weight 1 less belongs to the subgraph. |
---|
1889 | { |
---|
1890 | store=rtm_fix[j]-rmt_fix[j]-rtm_fix[j-1]; |
---|
1891 | } |
---|
1892 | // It is connected to two vertices V (which has weight one less) and W. |
---|
1893 | if (store>0) |
---|
1894 | { |
---|
1895 | if (j==stricttransforms[k]+1) // V is p[k] (to which the subgraph was glued) and W is q[k], the |
---|
1896 | { // vertex of weight one less, to which p[k] is connected. |
---|
1897 | rtm_var[j]=rtm_var[p[k]]+rtm_var[q[k]]; // In this case the subgraph separates p[k] and q[k]! |
---|
1898 | } |
---|
1899 | if (j==stricttransforms[k]+2) // V belongs to the subgraph (it is the vertex considerd in the |
---|
1900 | { // previous case!), and W is p[k] or q[k]. In this case the subgraph separates p[k] and q[k]. |
---|
1901 | if (store==rtm_fix[p[k]]) // If W=p[k], ... |
---|
1902 | { |
---|
1903 | rtm_var[j]=rtm_var[j-1]+rtm_var[p[k]]; |
---|
1904 | } |
---|
1905 | else // else W=q[k] ... . |
---|
1906 | { |
---|
1907 | rtm_var[j]=rtm_var[j-1]+rtm_var[q[k]]; // separates p[k] and q[k]. |
---|
1908 | } |
---|
1909 | } |
---|
1910 | if (j>stricttransforms[k]+2) // V belongs to the subgraph and W either does as well or is p[k]. |
---|
1911 | { |
---|
1912 | l=j-2; |
---|
1913 | while (store<rtm_fix[l]) // Find the second vertex W which is connected to which it is |
---|
1914 | { // connected. It has total multipl. = store! |
---|
1915 | if (l>stricttransforms[k]+1) // If l-1 belongs to the graph, then reduce l. |
---|
1916 | { |
---|
1917 | l=l-1; |
---|
1918 | } |
---|
1919 | else // If l-1 is stricttransform[k], hence isn't in the graph, then reducing l gives p[k]. |
---|
1920 | { // If l=p[k] and still store<rtm_fix[l], then j must be connected to q[k]! |
---|
1921 | if (l==stricttransforms[k]+1) |
---|
1922 | { |
---|
1923 | l=p[k]; |
---|
1924 | } |
---|
1925 | else |
---|
1926 | { |
---|
1927 | l=q[k]; |
---|
1928 | } |
---|
1929 | } |
---|
1930 | } |
---|
1931 | rtm_var[j]=rtm_var[j-1]+rtm_var[l]; |
---|
1932 | } |
---|
1933 | } |
---|
1934 | // It is only connected to one vertex V, which then must be the one of weight one less. |
---|
1935 | else |
---|
1936 | { |
---|
1937 | if (j==stricttransforms[k]+1) // V is p[k], the vertex, to which the subgraph was glued. |
---|
1938 | { |
---|
1939 | rtm_var[j]=rtm_var[p[k]]; |
---|
1940 | } |
---|
1941 | else |
---|
1942 | { |
---|
1943 | rtm_var[j]=rtm_var[j-1]; // V is belongs already to the subgraph. |
---|
1944 | } |
---|
1945 | } |
---|
1946 | } |
---|
1947 | } |
---|
1948 | return(rtm_var); |
---|
1949 | } |
---|
1950 | |
---|
1951 | |
---|
1952 | static proc find_connection_point (intvec v, int k) |
---|
1953 | "USAGE: find_connection_point(v,c); where v is an intvec, and k is an integer |
---|
1954 | RETURN: The largest integer i>k, such that v[i]=1, or 0 if no such i exists. |
---|
1955 | NOTE: This procedure is only for internal use; it is called in totalmultiplicities. |
---|
1956 | " |
---|
1957 | { |
---|
1958 | for (int i=size(v)-1;i>=k+1;i--) |
---|
1959 | { |
---|
1960 | if (v[i]==1) |
---|
1961 | { |
---|
1962 | return(i); |
---|
1963 | } |
---|
1964 | } |
---|
1965 | return(0); |
---|
1966 | } |
---|
1967 | |
---|
1968 | static proc find_connection_points (intmat resgr, int k) |
---|
1969 | "USAGE: find_connection_points(resgr,k); where resgr is an intmat, and k is an integer |
---|
1970 | RETURN: list of two intvecs ctps and ctpswts, where ctps contains all integers i!=k, such |
---|
1971 | that resgr[k,i]=1, and ctpswts contains for each such i the weight resgr[i,i]. |
---|
1972 | NOTE: This procedure is only for internal use; it is called in totalmultiplicities. |
---|
1973 | " |
---|
1974 | { |
---|
1975 | intvec ctps; |
---|
1976 | intvec ctpswts; |
---|
1977 | int j=1; |
---|
1978 | for (int i=1;i<=ncols(resgr);i++) |
---|
1979 | { |
---|
1980 | if ((resgr[k,i]==1) and (i!=k)) |
---|
1981 | { |
---|
1982 | ctps[j]=i; |
---|
1983 | ctpswts[j]=resgr[i,i]; |
---|
1984 | j++; |
---|
1985 | } |
---|
1986 | } |
---|
1987 | return(list(ctps,ctpswts)); |
---|
1988 | } |
---|
1989 | |
---|
1990 | static proc find_lower_connection_points (intmat resgr, int k) |
---|
1991 | "USAGE: find_lower_connection_points(resgr,k); where resgr is an intmat, and k is an integer |
---|
1992 | ASSUME: resgr is the resolutiongraph of an IRREDUCIBLE curve singularity and k<ncols(resgr). |
---|
1993 | RETURN: intvec, which contains the weights of points of lower weight than k, to which |
---|
1994 | the point of weight k in resgr is connected, and 0 if no such point exists. |
---|
1995 | NOTE: This procedure is only for internal use; it is called in totalmultiplicities. |
---|
1996 | " |
---|
1997 | { |
---|
1998 | intvec ctps=find_connection_points(resgr,k)[2]; |
---|
1999 | intvec lower_ctps; |
---|
2000 | int i=1; |
---|
2001 | while ((ctps[i]<k) and (ctps[i]>0)) |
---|
2002 | { |
---|
2003 | lower_ctps[i]=ctps[i]; |
---|
2004 | i++; |
---|
2005 | } |
---|
2006 | return(lower_ctps); |
---|
2007 | } |
---|
2008 | |
---|
2009 | |
---|
2010 | static proc euclides(int a,int b) // Procedure of Fernando. |
---|
2011 | "USAGE: euclides(m,n);where m,n are integers. |
---|
2012 | RETURN: The divisor, the quotients and the remainders of the euclidean algorithm performing for m and n. |
---|
2013 | NOTE: This procedure is only for internal use; it is called in multseq2charexp. |
---|
2014 | " |
---|
2015 | { |
---|
2016 | int c=a div b;//we compute in c the integer division between a and b. |
---|
2017 | int r=a mod b;//in r the remainer of the division between a and b |
---|
2018 | intvec dividendo=c;// we define the intvec of the dividens and we initialize it to c |
---|
2019 | intvec remain=r;// we define the intvec of the remainders and we initialize it to r |
---|
2020 | a=b;//change a to b |
---|
2021 | b=r;// and b to r |
---|
2022 | |
---|
2023 | while(r<>0)// while the current remainder is diferent to 0 we make the same af before |
---|
2024 | { |
---|
2025 | c=a div b; |
---|
2026 | r=a mod b; |
---|
2027 | dividendo =dividendo,c; |
---|
2028 | if(r<>0) |
---|
2029 | { |
---|
2030 | remain=remain,r; |
---|
2031 | } |
---|
2032 | a=b; |
---|
2033 | b=r; |
---|
2034 | } |
---|
2035 | list L=dividendo,remain;//we put in a list all the dividens and all the remainders |
---|
2036 | return(L);// and return the list |
---|
2037 | } |
---|
2038 | |
---|
2039 | |
---|
2040 | |
---|
2041 | static proc tau_es_qh (poly f) |
---|
2042 | "USAGE: tau_es_qh(f), poly f |
---|
2043 | RETURN: int, equisingular Tjurina number |
---|
2044 | NOTE: This procedure is only for internal use; it is called in Tau_es. |
---|
2045 | " |
---|
2046 | { |
---|
2047 | intvec qh=qhweight(f); |
---|
2048 | if (qh[1]==0) |
---|
2049 | { |
---|
2050 | ERROR("Input is not quasi-homogenous."); |
---|
2051 | } |
---|
2052 | else |
---|
2053 | { |
---|
2054 | int d_f = deg(f,qh); |
---|
2055 | list Tj=Tjurina(f,1); |
---|
2056 | ideal tj=Tj[2]; |
---|
2057 | int Taues=size(tj); |
---|
2058 | for (int i=1;i<=size(tj);i++) |
---|
2059 | { |
---|
2060 | if (deg(tj[i],qh)>=d_f) |
---|
2061 | { |
---|
2062 | Taues--; |
---|
2063 | } |
---|
2064 | } |
---|
2065 | } |
---|
2066 | return(Taues); |
---|
2067 | } |
---|
2068 | |
---|
2069 | |
---|
2070 | static proc move_row (intmat M, int i,j) |
---|
2071 | "USAGE: move_row(M,i,j), intmat M, int i,j |
---|
2072 | RETURN: intmat, the matrix M with j-th row now i-th row and the remaining rows moved accordingly. |
---|
2073 | NOTE: This procedure is only for internal use; it is called in move_row_col. |
---|
2074 | " |
---|
2075 | { |
---|
2076 | if ((i>nrows(M)) or (j>nrows(M))) |
---|
2077 | { |
---|
2078 | ERROR("The matrix has not enough rows."); |
---|
2079 | } |
---|
2080 | if (i==j) |
---|
2081 | { |
---|
2082 | return(M); |
---|
2083 | } |
---|
2084 | if (i>1) |
---|
2085 | { |
---|
2086 | intmat N[nrows(M)+1][ncols(M)]=intvec(M[1..i-1,1..ncols(M)]),intvec(M[j,1..ncols(M)]),intvec(M[i..nrows(M),1..ncols(M)]); |
---|
2087 | } |
---|
2088 | if (i==1) |
---|
2089 | { |
---|
2090 | intmat N[nrows(M)+1][ncols(M)]=intvec(M[j,1..ncols(M)]),intvec(M[i..nrows(M),1..ncols(M)]); |
---|
2091 | } |
---|
2092 | if (i<j) |
---|
2093 | { |
---|
2094 | N=delete_row(N,j+1); |
---|
2095 | } |
---|
2096 | if (i>j) |
---|
2097 | { |
---|
2098 | N=delete_row(N,j); |
---|
2099 | } |
---|
2100 | return(N); |
---|
2101 | } |
---|
2102 | |
---|
2103 | static proc move_col (intmat M, int i,j) |
---|
2104 | "USAGE: move_col(M,i,j), intmat M, int i,j |
---|
2105 | RETURN: intmat, the matrix M with j-th column now i-th column and the remaining columns moved accordingly. |
---|
2106 | NOTE: This procedure is only for internal use; it is called in move_row_col. |
---|
2107 | " |
---|
2108 | { |
---|
2109 | return(transpose(move_row(transpose(M),i,j))); |
---|
2110 | } |
---|
2111 | |
---|
2112 | static proc move_row_col (intmat M,int i,j) |
---|
2113 | "USAGE: move_row(M,i,j), intmat M, int i,j |
---|
2114 | RETURN: intmat, the matrix M with j-th row/column now i-th row/column and the remaining |
---|
2115 | rows and columns moved accordingly. |
---|
2116 | NOTE: This procedure is only for internal use; it is called in totalmultiplicities. |
---|
2117 | " |
---|
2118 | { |
---|
2119 | return(move_col(move_row(M,i,j),i,j)); |
---|
2120 | } |
---|
2121 | |
---|
2122 | |
---|
2123 | static proc delete_row (intmat M,int i) |
---|
2124 | "USAGE: delete_row(M,i); M intmat, i int |
---|
2125 | RETURN: intmat, which is derived from M by removing the ith row |
---|
2126 | NOTE: This procedure is only for internal use; it is called in move_row and tau_es2. |
---|
2127 | " |
---|
2128 | { |
---|
2129 | if ((i!=1) and (i!=nrows(M))) |
---|
2130 | { |
---|
2131 | return(intmat(intvec(M[1..i-1,1..ncols(M)],M[i+1..nrows(M),1..ncols(M)]),nrows(M)-1,ncols(M))); |
---|
2132 | } |
---|
2133 | if (i==1) |
---|
2134 | { |
---|
2135 | return(intmat(intvec(M[2..nrows(M),1..ncols(M)]),nrows(M)-1,ncols(M))); |
---|
2136 | } |
---|
2137 | else |
---|
2138 | { |
---|
2139 | return(intmat(intvec(M[1..nrows(M)-1,1..ncols(M)]),nrows(M)-1,ncols(M))); |
---|
2140 | } |
---|
2141 | } |
---|
2142 | |
---|
2143 | static proc delete_col (intmat M, int i) |
---|
2144 | "USAGE: delete_col(M,i); M intmat, i int |
---|
2145 | RETURN: intmat, which is derived from M by removing the ith column |
---|
2146 | NOTE: This procedure is only for internal use; it is called in tau_es. |
---|
2147 | " |
---|
2148 | { |
---|
2149 | return(transpose(delete_row(transpose(M),i))); |
---|
2150 | } |
---|
2151 | |
---|
2152 | |
---|
2153 | |
---|
2154 | static proc sort_branches (intmat contact, list graphs, list totmult, list multpl, int k,l) |
---|
2155 | "USAGE: sort_branches(K,L,M,N,k,l); K intmat, L,M,N lists, k,l integers |
---|
2156 | ASSUME: K = contact matrix of the branches of a curve, L = their resolutiongraphs, |
---|
2157 | M = their totalmultiplicities, N = their multiplicities |
---|
2158 | RETURN: list LL, LL[1] transformed K, LL[2..4] transformed L,M,N. |
---|
2159 | The procedure sorts the branches from k+1 to l according to descending contact with |
---|
2160 | with the branch k. |
---|
2161 | NOTE: This procedure is only for internal use; it is called in totalmultiplicities. |
---|
2162 | " |
---|
2163 | { |
---|
2164 | intvec max; |
---|
2165 | for (int i=k+1;i<=l;i++) |
---|
2166 | { |
---|
2167 | // Find the last graph max between i and l which has maximal contact with k. |
---|
2168 | max=max_in_intvec(intvec(contact[k,i..l]),1); |
---|
2169 | max[2]=max[2]+i-1; |
---|
2170 | if (i!=max[2]) // If max is not i, then move max to position i! |
---|
2171 | { |
---|
2172 | graphs=insert(graphs,graphs[max[2]],i-1); |
---|
2173 | graphs=delete(graphs,max[2]+1); |
---|
2174 | totmult=insert(totmult,totmult[max[2]],i-1); |
---|
2175 | totmult=delete(totmult,max[2]+1); |
---|
2176 | multpl=insert(multpl,multpl[max[2]],i-1); |
---|
2177 | multpl=delete(multpl,max[2]+1); |
---|
2178 | contact=move_row_col(contact,i,max[2]); |
---|
2179 | } |
---|
2180 | } |
---|
2181 | return(list(contact,graphs,totmult,multpl)); |
---|
2182 | } |
---|
2183 | |
---|
2184 | |
---|
2185 | static proc find_last_non_one (intvec v,int k) |
---|
2186 | "USAGE: find_last_non_one (v,k); intvec v, int k |
---|
2187 | RETURN: int i, the last index i>=k such that v[i]!=1, or k-1 if no such i exists. |
---|
2188 | NOTE: This procedure is only for internal use; it is called in tau_es2. |
---|
2189 | " |
---|
2190 | { |
---|
2191 | int i=size(v); |
---|
2192 | while (i>=k) |
---|
2193 | { |
---|
2194 | if (v[i]!=1) |
---|
2195 | { |
---|
2196 | return(i); |
---|
2197 | } |
---|
2198 | else |
---|
2199 | { |
---|
2200 | i--; |
---|
2201 | } |
---|
2202 | } |
---|
2203 | return(i); |
---|
2204 | } |
---|
2205 | |
---|
2206 | static proc intmat_minus_one (intmat M) |
---|
2207 | "USAGE: intmat_minus_one(M); intmat M |
---|
2208 | RETURN: intmat, all non-zero entries of M deminuished by 1. |
---|
2209 | NOTE: This procedure is only for internal use; it is called in tau_es2. |
---|
2210 | " |
---|
2211 | { |
---|
2212 | int i,j; |
---|
2213 | for (i=1;i<=nrows(M);i++) |
---|
2214 | { |
---|
2215 | for (j=1;j<=ncols(M);j++) |
---|
2216 | { |
---|
2217 | if (M[i,j]!=0) |
---|
2218 | { |
---|
2219 | M[i,j]=M[i,j]-1; |
---|
2220 | } |
---|
2221 | } |
---|
2222 | } |
---|
2223 | return(M); |
---|
2224 | } |
---|
2225 | |
---|
2226 | proc proximitymatrix (def INPUT) |
---|
2227 | "USAGE: proximitymatrix(INPUT); INPUT poly or list or intmat |
---|
2228 | ASSUME: INPUT is either a REDUCED bivariate polynomial defining a plane curve singularity, |
---|
2229 | or the output of @code{hnexpansion(f[,\"ess\"])}, or the list @code{hne} in |
---|
2230 | the ring created by @code{hnexpansion(f[,\"ess\"])}, or the output of |
---|
2231 | @code{develop(f)} resp. of @code{extdevelop(f,n)}, or a list containing |
---|
2232 | the contact matrix and a list of integer vectors with the characteristic exponents |
---|
2233 | of the branches of a plane curve singularity, or an integer vector containing the |
---|
2234 | characteristic exponents of an irreducible plane curve singularity, or the resolution |
---|
2235 | graph of a plane curve singularity (i.e. the output of resolutiongraph or |
---|
2236 | the first entry in the output of totalmultiplicities). |
---|
2237 | RETURN: list, of three integer matrices. The first one is the proximity matrix of |
---|
2238 | the plane curve defined by the INPUT, i.e. the entry i,j is 1 if the |
---|
2239 | infinitely near point corresponding to row i is proximate to the infinitely |
---|
2240 | near point corresponding to row j. The second integer matrix is the incidence matrix of the |
---|
2241 | resolution graph of the plane curve. The entry on the diagonal in row i is -s-1 |
---|
2242 | if s is the number of points proximate to the infinitely near point corresponding |
---|
2243 | to the ith row in the matrix. The third integer matrix is the incidence matrix of |
---|
2244 | the Enriques diagram of the plane curve singularity, i.e. each row corresponds to |
---|
2245 | an infinitely near point in the minimal standard resolution, including the |
---|
2246 | strict transforms of the branches, the diagonal element gives |
---|
2247 | the level of the point, and the entry i,j is -1 if row i is proximate to row j. |
---|
2248 | NOTE: In case the Hamburger-Noether expansion of the curve f is needed |
---|
2249 | for other purposes as well it is better to calculate this first |
---|
2250 | with the aid of @code{hnexpansion} and use it as input instead of |
---|
2251 | the polynomial itself. |
---|
2252 | @* |
---|
2253 | If you are not sure whether the INPUT polynomial is reduced or not, use |
---|
2254 | @code{squarefree(INPUT)} as input instead. |
---|
2255 | @* |
---|
2256 | If the input is a smooth curve, then the output will consist of |
---|
2257 | three one-by-one zero matrices. |
---|
2258 | @* |
---|
2259 | For the definitions of the computed objects see e.g. the book |
---|
2260 | Eduardo Casas-Alvero, Singularities of Plane Curves. |
---|
2261 | SEE ALSO: develop, hnexpansion, totalmultiplicities, alexanderpolynomial |
---|
2262 | EXAMPLE: example proximitymatrix; shows an example |
---|
2263 | " |
---|
2264 | { |
---|
2265 | ///////// Decide on the type of input. ////////// |
---|
2266 | if (typeof(INPUT)=="intmat") |
---|
2267 | { |
---|
2268 | intmat resgr=INPUT; |
---|
2269 | } |
---|
2270 | else |
---|
2271 | { |
---|
2272 | intmat resgr=totalmultiplicities(INPUT)[1]; |
---|
2273 | } |
---|
2274 | //////// Deal with the case of a smooth curve //////////////// |
---|
2275 | if (size(resgr)==1) |
---|
2276 | { |
---|
2277 | return(list(intmat(intvec(1),1,1),intmat(intvec(-1),1,1),intmat(intvec(0),1,1))); |
---|
2278 | } |
---|
2279 | //////// Calculate the proximity resolution graph //////////// |
---|
2280 | int i,j; |
---|
2281 | int n=nrows(resgr); |
---|
2282 | intvec diag; // Diagonal of the Enriques diagram. |
---|
2283 | int si; // number of points proximate to the point corresponding to the ith row |
---|
2284 | // Adjust the weights of the nodes corresponding to strict transforms so |
---|
2285 | // as if there had been one additional blow up. |
---|
2286 | for (i=1;i<=n;i++) |
---|
2287 | { |
---|
2288 | // Check if the row corresponds to an exceptional divisor ... |
---|
2289 | if (resgr[i,i]<0) |
---|
2290 | { |
---|
2291 | j=1; |
---|
2292 | while ((resgr[i,j]==0) or (i==j)) |
---|
2293 | { |
---|
2294 | j++; |
---|
2295 | } |
---|
2296 | resgr[i,i]=resgr[j,j]+1; |
---|
2297 | } |
---|
2298 | } |
---|
2299 | // Set the weights in the resolution graph to the blowing up level in the resolution. |
---|
2300 | for (i=1;i<=n;i++) |
---|
2301 | { |
---|
2302 | resgr[i,i]=resgr[i,i]-1; |
---|
2303 | diag[i]=resgr[i,i]; // The level of the corresponding infinitely near point. |
---|
2304 | } |
---|
2305 | // Replace the weights in the resolution graph by |
---|
2306 | // -s-1, where s is the number of points which are proximate to the point. |
---|
2307 | for (i=1;i<=n;i++) |
---|
2308 | { |
---|
2309 | si=-1; |
---|
2310 | // Find the points of higher weight which are connected to the ith row. |
---|
2311 | for (j=i+1;j<=n;j++) |
---|
2312 | { |
---|
2313 | // If the point in row j is connected to the ith row, then all the points of |
---|
2314 | // weight resgr[i,i]+1 to weight resgr[j,j] in the corresponding subgraph |
---|
2315 | // are proximate to the point of the ith row. This number is just resgr[j,j]-resgr[i,i]. |
---|
2316 | if ((resgr[i,j]!=0) and (resgr[j,j]>0)) |
---|
2317 | { |
---|
2318 | si=si-(resgr[j,j]-resgr[i,i]); |
---|
2319 | } |
---|
2320 | } |
---|
2321 | resgr[i,i]=si; |
---|
2322 | } |
---|
2323 | /////////////// Calculate the proximity matrix /////////////////// |
---|
2324 | intmat proximity=proxgauss(resgr); |
---|
2325 | intmat enriques=proximity; |
---|
2326 | for (i=1;i<=nrows(enriques);i++) |
---|
2327 | { |
---|
2328 | enriques[i,i]=diag[i]; |
---|
2329 | } |
---|
2330 | return(list(proximity,resgr,enriques)); |
---|
2331 | } |
---|
2332 | example |
---|
2333 | { |
---|
2334 | "EXAMPLE:"; |
---|
2335 | echo=2; |
---|
2336 | ring r=0,(x,y),ls; |
---|
2337 | poly f1=(y2-x3)^2-4x5y-x7; |
---|
2338 | poly f2=y2-x3; |
---|
2339 | poly f3=y3-x2; |
---|
2340 | list proximity=proximitymatrix(f1*f2*f3); |
---|
2341 | /// The proximity matrix P /// |
---|
2342 | print(proximity[1]); |
---|
2343 | /// The proximity resolution graph N /// |
---|
2344 | print(proximity[2]); |
---|
2345 | /// They satisfy N=-transpose(P)*P /// |
---|
2346 | print(-transpose(proximity[1])*proximity[1]); |
---|
2347 | /// The incidence matrix of the Enriques diagram /// |
---|
2348 | print(proximity[3]); |
---|
2349 | /// If M is the matrix of multiplicities and TM the matrix of total |
---|
2350 | /// multiplicities of the singularity, then M=P*TM. |
---|
2351 | /// We therefore calculate the (total) multiplicities. Note that |
---|
2352 | /// they have to be slightly extended. |
---|
2353 | list MULT=extend_multiplicities(totalmultiplicities(f1*f2*f3)); |
---|
2354 | intmat TM=MULT[1]; // Total multiplicites. |
---|
2355 | intmat M=MULT[2]; // Multiplicities. |
---|
2356 | /// Check: M-P*TM=0. |
---|
2357 | M-proximity[1]*TM; |
---|
2358 | /// Check: inverse(P)*M-TM=0. |
---|
2359 | intmat_inverse(proximity[1])*M-TM; |
---|
2360 | } |
---|
2361 | |
---|
2362 | static proc addmultiplrows (intmat M, int i, int j, int ki, int kj) |
---|
2363 | "USAGE: addmultiplrows(M,i,j,ki,kj); intmat M, int i,j,ki,kj |
---|
2364 | RETURN: intmat, replaces the j-th row in M by ki-times the i-th row plus |
---|
2365 | kj times the j-th |
---|
2366 | NOTE: This procedure is only for internal use; it is called in intmat_inverse |
---|
2367 | and proxgauss. |
---|
2368 | " |
---|
2369 | { |
---|
2370 | int k=ncols(M); |
---|
2371 | M[j,1..k]=kj*M[j,1..k]+ki*M[i,1..k]; |
---|
2372 | return(M); |
---|
2373 | } |
---|
2374 | |
---|
2375 | |
---|
2376 | static proc proxgauss (intmat M) |
---|
2377 | "USAGE: proxgauss(M); intmat M |
---|
2378 | ASSUME: M is the output of proximity_resgr |
---|
2379 | RETURN: intmat, replaces the j-th row in M by ki-times the i-th row plus |
---|
2380 | kj times the j-th |
---|
2381 | NOTE: This procedure is only for internal use; it is called in intmat_inverse. |
---|
2382 | " |
---|
2383 | { |
---|
2384 | int i; |
---|
2385 | int n=nrows(M); |
---|
2386 | if (n==1) |
---|
2387 | { |
---|
2388 | M[1,1]=1; |
---|
2389 | return(M); |
---|
2390 | } |
---|
2391 | else |
---|
2392 | { |
---|
2393 | if (M[n,n]<0) |
---|
2394 | { |
---|
2395 | M=addmultiplrows(M,n,n,-1,0); |
---|
2396 | } |
---|
2397 | for (i=n-1;i>=1;i--) |
---|
2398 | { |
---|
2399 | if (M[i,n]!=0) |
---|
2400 | { |
---|
2401 | M=addmultiplrows(M,n,i,-M[i,n],M[n,n]); |
---|
2402 | } |
---|
2403 | } |
---|
2404 | intmat N[n-1][n-1]=M[1..n-1,1..n-1]; |
---|
2405 | N=proxgauss(N); |
---|
2406 | M[1..n-1,1..n-1]=N[1..n-1,1..n-1]; |
---|
2407 | return(M); |
---|
2408 | } |
---|
2409 | } |
---|
2410 | |
---|
2411 | proc extend_multiplicities (list rg) |
---|
2412 | "USAGE: extend_multiplicities(rg); list rg |
---|
2413 | ASSUME: rg is the output of the procedure totalmultiplicities |
---|
2414 | RETURN: list, the first entry is an integer matrix containing the total |
---|
2415 | multiplicities and the second entry is an integer matrix containing |
---|
2416 | the multiplicies of the resolution given by rg, where the zeros |
---|
2417 | corresponding to the strict transforms of the branches of the curve |
---|
2418 | have been replaced by the (total) multiplicities of the infinitely near |
---|
2419 | point corresponding to one further blow up for each branch. |
---|
2420 | KEYWORDS: total multiplicities; multiplicity sequence; resolution graph |
---|
2421 | EXAMPLE: example extend_multiplicities; shows an example |
---|
2422 | " |
---|
2423 | { |
---|
2424 | intmat resgr,tm,mt=rg[1],rg[2],rg[3]; |
---|
2425 | int i,j; |
---|
2426 | for (i=1;i<=nrows(resgr);i++) |
---|
2427 | { |
---|
2428 | if (resgr[i,i]<0) |
---|
2429 | { |
---|
2430 | j=1; |
---|
2431 | while ((resgr[i,j]==0) or (i==j)) |
---|
2432 | { |
---|
2433 | j++; |
---|
2434 | } |
---|
2435 | tm[i,1..ncols(tm)]=tm[j,1..ncols(tm)]; |
---|
2436 | tm[i,-resgr[i,i]]=tm[i,-resgr[i,i]]+1; |
---|
2437 | mt[i,-resgr[i,i]]=1; |
---|
2438 | } |
---|
2439 | } |
---|
2440 | return(list(tm,mt)); |
---|
2441 | } |
---|
2442 | example |
---|
2443 | { |
---|
2444 | "EXAMPLE:"; |
---|
2445 | echo=2; |
---|
2446 | ring r=0,(x,y),ls; |
---|
2447 | poly f1=(y2-x3)^2-4x5y-x7; |
---|
2448 | poly f2=y2-x3; |
---|
2449 | poly f3=y3-x2; |
---|
2450 | // Calculate the resolution graph and the (total) multiplicities of f1*f2*f3. |
---|
2451 | list RESGR=totalmultiplicities(f1*f2*f3); |
---|
2452 | // Extend the (total) multiplicities. |
---|
2453 | list MULT=extend_multiplicities(RESGR); |
---|
2454 | // Compare the total multiplicities. |
---|
2455 | RESGR[2]; |
---|
2456 | MULT[1]; |
---|
2457 | // Compare the multiplicities. |
---|
2458 | RESGR[3]; |
---|
2459 | MULT[2]; |
---|
2460 | } |
---|
2461 | |
---|
2462 | proc intmat_inverse (intmat M) |
---|
2463 | "USAGE: intmat_inverse(M); intmat M |
---|
2464 | ASSUME: M is a lower triangular integer matrix with diagonal entries 1 or -1 |
---|
2465 | RETURN: intmat, the inverse of M |
---|
2466 | KEYWORDS: integer matrix, inverse |
---|
2467 | EXAMPLE: example intmat_inverse; shows an example |
---|
2468 | " |
---|
2469 | { |
---|
2470 | int i,j,k; |
---|
2471 | int n=nrows(M); |
---|
2472 | intmat U[n][n]; |
---|
2473 | U=U+1; |
---|
2474 | for (i=1;i<=n;i++) |
---|
2475 | { |
---|
2476 | if (M[i,i]==-1) |
---|
2477 | { |
---|
2478 | M=addmultiplrows(M,i,i,-1,0); |
---|
2479 | U=addmultiplrows(U,i,i,-1,0); |
---|
2480 | } |
---|
2481 | for (j=i+1;j<=n;j++) |
---|
2482 | { |
---|
2483 | U=addmultiplrows(U,i,j,-M[j,i],M[i,i]); |
---|
2484 | M=addmultiplrows(M,i,j,-M[j,i],M[i,i]); |
---|
2485 | } |
---|
2486 | } |
---|
2487 | return(U); |
---|
2488 | } |
---|
2489 | example |
---|
2490 | { |
---|
2491 | "EXAMPLE:";echo=2; |
---|
2492 | intmat M[5][5]=1,0,0,0,0,1,1,0,0,0,2,1,1,0,0,3,1,1,1,0,4,1,1,1,1 ; |
---|
2493 | intmat U=intmat_inverse(M); |
---|
2494 | print(U); |
---|
2495 | print(U*M); |
---|
2496 | } |
---|