1 // ----------------------------------------------------------------------------
2 //  semic.cc
3 //  begin of file
4 //  Stephan Endrass, endrass@mathematik.uni-mainz.de
5 //  23.7.99
6 // ----------------------------------------------------------------------------
7 
8 #define SEMIC_CC
9 
10 
11 
12 
13 #include "kernel/mod2.h"
14 
15 #ifdef HAVE_SPECTRUM
16 
17 #ifdef  SEMIC_PRINT
18 #ifndef SEMIC_IOSTREAM
19 #include <stdio.h>
20 #else
21 #include <iostream.h>
22 #endif
23 #endif
24 
25 #include <string.h>
26 
27 #include "misc/intvec.h"
28 #include "misc/mylimits.h"
29 #include "kernel/spectrum/GMPrat.h"
30 #include "kernel/spectrum/semic.h"
31 
32 // ----------------------------------------------------------------------------
33 //  Copy constructor for  spectrum
34 // ----------------------------------------------------------------------------
35 
spectrum(const spectrum & spec)36 spectrum::spectrum( const spectrum &spec )
37 {
38     copy_deep( spec );
39 }
40 
41 // ----------------------------------------------------------------------------
42 //  Destructor for  spectrum
43 // ----------------------------------------------------------------------------
44 
~spectrum()45 spectrum::~spectrum( )
46 {
47     copy_delete( );
48 }
49 
50 // ----------------------------------------------------------------------------
51 //  Allocate memory for a spectrum of  k  numbers
52 // ----------------------------------------------------------------------------
53 
copy_new(int k)54 void spectrum::copy_new( int k )
55 {
56     if( k > 0 )
57     {
58         s = new Rational[k];
59         w = new int[k];
60 
61         #ifndef SING_NDEBUG
62         if( s == (Rational*)NULL || w == (int*)NULL )
63         {
64             #ifdef SEMIC_PRINT
65             #ifdef SEMIC_IOSTREAM
66                 cerr << "spectrum::copy_new(" << k << ")" << endl;
67                 cerr << "    returned ZERO!!!" << endl;
68                 cerr << "    exit..." << endl;
69             #else
70                 fprintf( stderr,"spectrum::copy_new( %d )\n",k );
71                 fprintf( stderr,"    returned ZERO!!!\n" );
72                 fprintf( stderr,"    exit...\n" );
73             #endif
74             #endif
75         }
76         #endif
77     }
78     else if( k == 0 )
79     {
80         s = (Rational*)NULL;
81         w = (int*)NULL;
82     }
83     else if( k < 0 )
84     {
85         #ifdef SEMIC_PRINT
86         #ifdef SEMIC_IOSTREAM
87                 cerr << "spectrum::copy_new(" << k << ")";
88                 cerr << ": k < 0 ..." << endl;
89         #else
90                 fprintf( stderr,"spectrum::copy_new( %d )",k );
91                 fprintf( stderr,": k < 0 ...\n" );
92         #endif
93         #endif
94 
95         exit( 1 );
96     }
97 }
98 
99 // ----------------------------------------------------------------------------
100 //  Initialize a  spectrum  deep from another  spectrum
101 // ----------------------------------------------------------------------------
102 
copy_deep(const spectrum & spec)103 void spectrum::copy_deep( const spectrum &spec )
104 {
105     mu = spec.mu;
106     pg = spec.pg;
107     n  = spec.n;
108 
109     copy_new( n );
110 
111     for( int i=0; i<n; i++ )
112     {
113         s[i] = spec.s[i];
114         w[i] = spec.w[i];
115     }
116 }
117 
118 // ----------------------------------------------------------------------------
119 //  operator  =  for  spectrum
120 // ----------------------------------------------------------------------------
121 
operator =(const spectrum & spec)122 spectrum spectrum::operator = ( const spectrum &spec )
123 {
124     copy_delete( );
125     copy_deep( spec );
126 
127     return *this;
128 }
129 
130 // ----------------------------------------------------------------------------
131 //  add the two spectra  s1  and  s2  and return their sum
132 // ----------------------------------------------------------------------------
133 
operator +(const spectrum & s1,const spectrum & s2)134 spectrum  operator + ( const spectrum &s1,const spectrum &s2 )
135 {
136     int i1=0, i2=0, i3=0;
137 
138     spectrum result;
139 
140     do
141     {
142         if( i1 >= s1.n )
143         {
144             i2++;
145         }
146         else if( i2 >= s2.n )
147         {
148             i1++;
149         }
150         else if( s1.s[i1] < s2.s[i2] )
151         {
152             i1++;
153         }
154         else if( s1.s[i1] == s2.s[i2] )
155         {
156             i1++;
157             i2++;
158         }
159         else
160         {
161             i2++;
162         }
163         i3++;
164     }
165     while( i1 < s1.n || i2 < s2.n );
166 
167     result.copy_new( i3 );
168     result.n = i3;
169 
170     i1 = i2 = i3 = 0;
171 
172     do
173     {
174         if( i1 >= s1.n )
175         {
176             result.s[i3] = s2.s[i2];
177             result.w[i3] = s2.w[i2];
178             i2++;
179         }
180         else if( i2 >= s2.n )
181         {
182             result.s[i3] = s1.s[i1];
183             result.w[i3] = s1.w[i1];
184             i1++;
185         }
186         else if( s1.s[i1] < s2.s[i2] )
187         {
188             result.s[i3] = s1.s[i1];
189             result.w[i3] = s1.w[i1];
190             i1++;
191           }
192         else if( s1.s[i1] == s2.s[i2] )
193         {
194             result.s[i3] = s1.s[i1];
195             result.w[i3] = s1.w[i1] + s2.w[i2];
196             i1++;
197             i2++;
198         }
199         else
200         {
201             result.s[i3] = s2.s[i2];
202             result.w[i3] = s2.w[i2];
203             i2++;
204         }
205         i3++;
206     }
207     while( i1 < s1.n || i2 < s2.n );
208 
209     result.mu = s1.mu + s2.mu;
210     result.pg = s1.pg + s2.pg;
211 
212     return  result;
213 }
214 
215 // ----------------------------------------------------------------------------
216 //  multiply the multiplicities of the spectrum numbers of  a  with  m
217 // ----------------------------------------------------------------------------
218 
operator *(int k,const spectrum & spec)219 spectrum operator * ( int k,const spectrum &spec )
220 {
221     if( k == 0 )
222     {
223         spectrum result;
224 
225         return  result;
226     }
227     else
228     {
229         spectrum result( spec );
230 
231         result.mu *= k;
232         result.pg *= k;
233 
234         for( int i=0; i<result.n; i++ )
235         {
236             result.w[i] *= k;
237         }
238 
239         return  result;
240     }
241 }
242 
243 // ----------------------------------------------------------------------------
244 //  Print a  spectrum
245 // ----------------------------------------------------------------------------
246 
247 #ifdef SEMIC_PRINT
248 
operator <<(ostream & s,const spectrum & spec)249 ostream & operator << ( ostream &s,const spectrum &spec )
250 {
251     for( int i=0; i<spec.n; i++ )
252     {
253         if( i>0 )
254         {
255             #ifdef SEMIC_STDOUT
256                 s << "+";
257             #else
258                 fprintf( stdout,"+" );
259             #endif
260         }
261 
262         #ifdef SEMIC_STDOUT
263             s << spec.w[i] << "*t^";
264         #else
265             fprintf( stdout,"%d*t^",spec.w[i] );
266         #endif
267 
268         s << spec.s[i];
269     }
270 
271     return s;
272 }
273 #endif
274 
275 // ----------------------------------------------------------------------------
276 //  Add a subspectrum with multiplicity  k  (faster than '+')
277 // ----------------------------------------------------------------------------
278 
add_subspectrum(spectrum & a,int k)279 int    spectrum::add_subspectrum( spectrum &a,int k )
280 {
281     int i,j;
282     for( i=0, j=0; i<n; i++ )
283     {
284         if( s[i] == a.s[j] )
285         {
286             w[i] += k*a.w[j];
287             j++;
288         }
289     }
290 
291     return ( j == a.n ? TRUE : FALSE );
292 }
293 
294 // ----------------------------------------------------------------------------
295 //  set  *alpha  to the next spectrum number strictly bigger than  *alpha
296 //  returns: TRUE, if such a spectrum number exists
297 //           FALSE otherwise
298 // ----------------------------------------------------------------------------
299 
next_number(Rational * alpha)300 int    spectrum::next_number( Rational *alpha )
301 {
302     int i=0;
303 
304     while( i < n && *alpha >= s[i]  )
305     {
306         i++;
307     }
308 
309     if( i < n )
310     {
311         *alpha = s[i];
312         return TRUE;
313     }
314     else
315     {
316         return FALSE;
317     }
318 }
319 
320 // ----------------------------------------------------------------------------
321 //  find the next interval on the real line of same length as
322 //  [*alpha1,*alpha2]  having a spectrum number as interval border
323 // ----------------------------------------------------------------------------
324 
next_interval(Rational * alpha1,Rational * alpha2)325 int     spectrum::next_interval( Rational *alpha1,Rational *alpha2 )
326 {
327     Rational zero( 0,1 );
328     Rational a1 = *alpha1;
329     Rational a2 = *alpha2;
330     Rational d  = *alpha2 - *alpha1;
331 
332     int    e1 = this->next_number( &a1 );
333     int    e2 = this->next_number( &a2 );
334 
335     if( e1 || e2 )
336     {
337         Rational d1 = a1 - *alpha1;
338         Rational d2 = a2 - *alpha2;
339 
340         if( d1 < d2 || d2 == zero )
341         {
342             *alpha1 = a1;
343             *alpha2 = a1 + d;
344         }
345         else
346         {
347             *alpha1 = a2 - d;
348             *alpha2 = a2;
349         }
350         return  TRUE;
351     }
352     else
353     {
354         return  FALSE;
355     }
356 }
357 
358 // ----------------------------------------------------------------------------
359 //  compute the numver of spectrum numbers in the inverval  [*alpha1,*alpha2]
360 // ----------------------------------------------------------------------------
361 
numbers_in_interval(Rational & alpha1,Rational & alpha2,interval_status status)362 int     spectrum::numbers_in_interval( Rational &alpha1,
363                 Rational &alpha2,interval_status status )
364 {
365     int count = 0;
366 
367     for( int i=0; i<n; i++ )
368     {
369         if( ( ( status == OPEN   || status == LEFTOPEN  ) &&
370               s[i] >  alpha1 ) ||
371             ( ( status == CLOSED || status == RIGHTOPEN ) &&
372               s[i] >= alpha1 ) )
373         {
374               if( ( ( status == OPEN   || status == RIGHTOPEN  ) &&
375                   s[i] <  alpha2 ) ||
376                 ( ( status == CLOSED || status == LEFTOPEN ) &&
377                   s[i] <= alpha2 ) )
378             {
379                 count += w[i];
380             }
381             else
382             {
383                 break;
384             }
385         }
386     }
387 
388     return count;
389 }
390 
391 // ----------------------------------------------------------------------------
392 //  find the maximal integer  k  such that  k*t is semicontinous
393 //  for the spectrum
394 // ----------------------------------------------------------------------------
395 
mult_spectrum(spectrum & t)396 int     spectrum::mult_spectrum( spectrum &t )
397 {
398     spectrum u = *this + t;
399 
400     Rational alpha1 = -2;
401     Rational alpha2 = -1;
402 
403     int      mult=MAX_INT_VAL,nthis,nt;
404 
405     while( u.next_interval( &alpha1,&alpha2 ) )
406     {
407         nt    = t.numbers_in_interval( alpha1,alpha2,LEFTOPEN );
408         nthis = this->numbers_in_interval( alpha1,alpha2,LEFTOPEN );
409 
410         if( nt != 0 )
411         {
412             mult = (nthis/nt < mult ? nthis/nt: mult );
413         }
414 
415     }
416 
417     return  mult;
418 }
419 
420 // ----------------------------------------------------------------------------
421 //  find the maximal integer  k  such that  k*t is semicontinous
422 //  for the spectrum (in the homogeneous sense)
423 // ----------------------------------------------------------------------------
424 
mult_spectrumh(spectrum & t)425 int     spectrum::mult_spectrumh( spectrum &t )
426 {
427     spectrum u = *this + t;
428 
429     Rational alpha1 = -2;
430     Rational alpha2 = -1;
431 
432     int      mult=MAX_INT_VAL,nthis,nt;
433 
434     while( u.next_interval( &alpha1,&alpha2 ) )
435     {
436         nt    = t.numbers_in_interval( alpha1,alpha2,LEFTOPEN );
437         nthis = this->numbers_in_interval( alpha1,alpha2,LEFTOPEN );
438 
439         if( nt != 0 )
440         {
441             mult = (nthis/nt < mult ? nthis/nt: mult );
442         }
443 
444         nt    = t.numbers_in_interval( alpha1,alpha2,OPEN );
445         nthis = this->numbers_in_interval( alpha1,alpha2,OPEN );
446 
447         if( nt != 0 )
448         {
449             mult = (nthis/nt < mult ? nthis/nt: mult );
450         }
451     }
452 
453     return  mult;
454 }
455 
456 // ----------------------------------------------------------------------------
457 //  Set the Milnor number
458 // ----------------------------------------------------------------------------
459 
460 /*
461 int spectrum::set_milnor( void )
462 {
463    mu = 0;
464 
465    for( int i=0; i<n; i++ )
466    {
467       mu += w[i];
468    }
469 
470    return  mu;
471 }
472 
473 // ----------------------------------------------------------------------------
474 //  Set the geometrical genus
475 // ----------------------------------------------------------------------------
476 
477 int spectrum::set_geometric_genus( void )
478 {
479    pg = 0;
480 
481    for( int i=0; i<n && s[i]<=1; i++ )
482    {
483       pg += w[i];
484    }
485    return  pg;
486 }
487 */
488 
489 #endif /* HAVE_SPECTRUM */
490 // ----------------------------------------------------------------------------
491 //  semic.cc
492 //  end of file
493 // ----------------------------------------------------------------------------
494