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