58 WORD EpfFind(PHEAD WORD *term, WORD *params)
61 WORD *t, *m, *r, n1 = 0, n2, min = -1, count, fac;
62 WORD *c1 = 0, *c2 = 0, sgn = 1;
64 UWORD *facto = (UWORD *)AT.WorkPointer;
67 if ( ( AT.WorkPointer = (WORD *)(facto + AM.MaxTal) ) > AT.WorkTop ) {
68 MLOCK(ErrorMessageLock);
70 MUNLOCK(ErrorMessageLock);
79 while ( *t != LEVICIVITA && t < tstop ) t += t[1];
80 if ( t >= tstop )
return(0);
82 while ( *m == LEVICIVITA && m < tstop ) { n1++; m += m[1]; }
84 if ( n1 <= (number+1) || n1 <= 1 )
return(0);
88 while ( m[1] == t[1] ) {
91 count = fac = n1 = n2 = t[1] - FUNHEAD;
100 if ( *m >= AM.OffsetIndex &&
101 ( ( *m >= AM.IndDum && AC.lDefDim == fac ) ||
103 indices[*m-AM.OffsetIndex].dimension == fac ) ) ) {
106 r++; m++; n1--; n2--;
110 if ( min < 0 || count < min ) {
112 c2 = m - fac - FUNHEAD;
115 if ( m >= mstop )
break;
118 }
while ( ( m = t + t[1] ) < mstop );
121 fac = type + FUNHEAD;
122 while ( *t != LEVICIVITA && t < tstop ) t += t[1];
123 while ( *t == LEVICIVITA && t < tstop && t[1] != fac ) t += t[1];
124 if ( t >= tstop )
return(0);
126 while ( *m == LEVICIVITA && m < tstop && m[1] == fac ) { n1++; m += m[1]; }
134 if ( min < 0 )
return(0);
137 fac = t[1] - FUNHEAD;
141 if ( number < 0 ) *m++ = 1;
143 n1 = n2 = t[1] - FUNHEAD;
149 if ( *m > *r ) { *c1++ = *r++; n2--; }
150 else if ( *m < *r ) { *c2++ = *m++; n1--; }
152 if ( *m < AM.OffsetIndex || ( *m < AM.IndDum &&
153 ( indices[*m-AM.OffsetIndex].dimension != fac ) ) ||
154 ( *m >= AM.IndDum && AC.lDefDim != fac ) ) {
155 *c1++ = *r++; *c2++ = *m++;
157 else {
if ( ( n1 ^ n2 ) & 1 ) sgn = -sgn; r++; m++; }
161 if ( n1 ) { NCOPY(c2,m,n1); }
162 else if ( n2 ) { NCOPY(c1,r,n2); }
165 while ( m < mstop ) *t++ = *m++;
167 while ( m < tstop ) *t++ = *m++;
168 *t++ = SUBEXPRESSION;
179 while ( m < mstop ) *t++ = *m++;
180 if ( Factorial(BHEAD fac,facto,&nfac) || Mully(BHEAD (UWORD *)tstop,&ncoef,facto,nfac) ) {
181 MLOCK(ErrorMessageLock);
183 MUNLOCK(ErrorMessageLock);
186 tstop += (ABS(ncoef))*2;
187 if ( sgn < 0 ) ncoef = -ncoef;
189 *tstop++ = (ncoef<0)?(ncoef-1):(ncoef+1);
190 *term = WORDDIF(tstop,term);
212 WORD EpfCon(PHEAD WORD *term, WORD *params, WORD num, WORD level)
215 WORD *kron, *perm, *termout, *tstop, size2;
216 WORD *m, *t, sizes, sgn = 0, i;
218 kron = AT.WorkPointer;
219 perm = (AT.WorkPointer += sizes);
220 termout = (AT.WorkPointer += sizes);
221 AT.WorkPointer = (WORD *)(((UBYTE *)(AT.WorkPointer)) + AM.MaxTer);
222 if ( AT.WorkPointer > AT.WorkTop ) {
223 MLOCK(ErrorMessageLock);
225 MUNLOCK(ErrorMessageLock);
229 if ( !(*params++) ) level--;
231 if ( !size2 )
goto DoOnce;
232 while ( ( sgn = EpfGen(size2,params,kron,perm,sgn) ) != 0 ) {
238 while ( t < tstop ) *m++ = *t++;
239 if ( t[2] != num || *t != SUBEXPRESSION ) {
240 MLOCK(ErrorMessageLock);
241 MesPrint(
"Serious error in EpfCon");
242 MUNLOCK(ErrorMessageLock);
251 if ( i ) { NCOPY(m,t,i); }
254 tstop = term + *term;
255 while ( t < tstop ) *m++ = *t++;
256 *termout = WORDDIF(m,termout);
258 if ( sgn < 0 ) *m = - *m;
262 AT.WorkPointer = termout + *termout;
263 if (
Generator(BHEAD termout,level) < 0 )
goto EpfCall;
266 AT.WorkPointer = kron;
269 if ( AM.tracebackflag ) {
270 MLOCK(ErrorMessageLock);
272 MUNLOCK(ErrorMessageLock);
282 WORD EpfGen(WORD number, WORD *inlist, WORD *kron, WORD *perm, WORD sgn)
286 in2 = inlist + number;
288 for ( i = 1; i < number; i += 2 ) {
294 if ( number <= 0 )
return(0);
299 while ( ( i -= 2 ) >= 0 ) {
300 if ( ( k = perm[i] ) != i ) {
306 if ( ( k = ( perm[i] += 2 ) ) < number ) {
311 for ( k = i + 2; k < number; k += 2 ) perm[k] = k;
331 WORD Trick(WORD *in,
TRACES *t)
336 switch ( t->eers[n] ) {
339 p = t->pepf + t->mepf[n];
352 p = t->pdel + t->mdel[n];
357 *in = *(t->pepf + t->mepf[n] + 2);
363 t->sign1 = - t->sign1;
364 *(t->pdel + t->mdel[n] + 1) = in[2];
369 t->sign1 = - t->sign1;
371 *(t->pdel + t->mdel[n]) = in[1];
375 in[2] = *(t->pdel + t->mdel[n] + 1);
383 return(--(t->eers[n]));
418 WORD Trace4no(WORD number, WORD *kron,
TRACES *t)
422 WORD retval, *stop, oldsign;
424 if ( ( number < 0 ) || ( number & 1 ) )
return(0);
426 if ( t->gamma5 == GAMMA5 )
return(0);
428 *kron++ = *t->inlist;
429 *kron++ = t->inlist[1];
435 WORD nhalf = number >> 1;
436 WORD ndouble = number * 2;
438 t->eers = p; p += nhalf;
439 t->mepf = p; p += nhalf;
440 t->mdel = p; p += nhalf;
441 t->pdel = p; p += number + nhalf;
442 t->pepf = p; p += ndouble;
443 t->e4 = p; p += number;
444 t->e3 = p; p += ndouble;
445 t->nt3 = p; p += nhalf;
446 t->nt4 = p; p += nhalf;
447 t->j3 = p; p += ndouble;
452 t->mdum = AM.mTraceDum;
461 t->stap = (t->step1)++;
463 t->eers[t->stap] = 0;
464 t->mepf[t->step1] = t->mepf[t->stap];
465 t->mdel[t->step1] = t->mdel[t->stap];
466 CallTrick:
while ( !Trick(t->inlist+t->kstep,t) ) {
468 t->step1 = (t->stap)--;
473 }
while ( t->kstep < (number-4) );
480 if ( ( t->gamma5 == GAMMA7 ) && ( t->gamm == -1 ) ) {
481 t->sign2 = - t->sign2;
483 else if ( ( t->gamma5 == GAMMA5 ) && ( t->gamm == 1 ) ) {
486 else if ( ( t->gamma5 == GAMMA1 ) && ( t->gamm == -1 ) ) {
489 p = t->pdel + t->mdel[t->step1];
490 *p++ = t->inlist[t->kstep+2];
491 *p++ = t->inlist[t->kstep+3];
505 ae = t->mepf[t->step1];
506 t->ad = t->mdel[t->step1]+2;
509 while ( ( ae -= 4 ) >= 0 ) {
510 if ( t->pepf[ae] > AM.mTraceDum && t->pepf[ae] <= t->mdum ) {
513 for ( i = 0; i < 3; i++ ) {
523 for ( i = 0; i < 4; i++ ) *p++ = *m++;
537 while ( t->a3 > 0 ) {
538 t->nt3[++(t->lc3)] = 0;
539 while ( ( t->nt3[t->lc3] = EpfGen(3,t->e3+t->a3-6,
540 t->pdel+t->ad,t->j3+6*t->lc3,oldsign = t->nt3[t->lc3]) ) == 0 ) {
541 if ( oldsign < 0 ) t->sign2 = - t->sign2;
543 NextE3:
if ( t->lc3 < 0 )
goto CallTrick;
548 if ( oldsign != t->nt3[t->lc3] ) t->sign2 = - t->sign2;
550 else if ( t->nt3[t->lc3] < 0 ) t->sign2 = - t->sign2;
557 while ( t->a4 > 4 ) {
558 t->nt4[++(t->lc4)] = 0;
559 while ( ( t->nt4[t->lc4] = EpfGen(4,t->e4+t->a4-8,
560 t->pdel+t->ad,t->j4+8*t->lc4,oldsign = t->nt4[t->lc4]) ) == 0 ) {
561 if ( oldsign < 0 ) t->sign2 = - t->sign2;
563 NextE4:
if ( t->lc4 < 0 )
goto NextE3;
568 if ( oldsign != t->nt4[t->lc4] ) t->sign2 = - t->sign2;
570 else if ( t->nt4[t->lc4] < 0 ) t->sign2 = - t->sign2;
584 *m++ = *p++; *m++ = *p++; *m++ = *p++; *m++ = *p++;
588 if ( t->sign2 < 0 ) retval = - retval;
590 for ( i = 0; i < t->ad; i++ ) *m++ = *p++;
596 stop = p + t->ad + 4;
598 while ( *p >= AM.mTraceDum && *p <= t->mdum ) {
607 else if ( m[1] == *p ) {
614 }
while ( m < stop );
618 else stop = p + t->ad;
619 while ( p < (stop-2) ) {
620 while ( *p >= AM.mTraceDum && *p <= t->mdum ) {
629 else if ( m[1] == *p ) {
636 }
while ( m < stop );
639 while ( *p >= AM.mTraceDum && *p <= t->mdum ) {
648 else if ( m[1] == *p ) {
655 }
while ( m < stop );
661 if ( number <= 2 )
return(0);
662 else {
goto NextE4; }
694 WORD Trace4(PHEAD WORD *term, WORD *params, WORD num, WORD level)
698 WORD *p, *m, number, i;
700 WORD j, minimum, minimum2, *min, *stopper;
701 OldW = AT.WorkPointer;
702 if ( AN.numtracesctack >= AN.intracestack ) {
703 number = AN.intracestack + 2;
704 t = (
TRACES *)Malloc1(number*
sizeof(
TRACES),
"TRACES-struct");
705 if ( AN.tracestack ) {
706 for ( i = 0; i < AN.intracestack; i++ ) { t[i] = AN.tracestack[i]; }
707 M_free(AN.tracestack,
"TRACES-struct");
710 AN.intracestack = number;
713 number = *params - 6;
714 if ( number < 0 || ( number & 1 ) || !params[5] )
return(0);
716 t = AN.tracestack + AN.numtracesctack;
719 t->finalstep = ( params[2] & 16 ) ? 1 : 0;
720 t->gamma5 = params[3];
721 if ( t->finalstep && t->gamma5 != GAMMA1 ) {
722 MLOCK(ErrorMessageLock);
723 MesPrint(
"Gamma5 not allowed in this option of the trace command");
724 MUNLOCK(ErrorMessageLock);
728 t->inlist = AT.WorkPointer;
729 t->accup = t->accu = t->inlist + number;
730 t->perm = t->accu + (number*2);
731 t->eers = t->perm + number;
732 if ( ( AT.WorkPointer += 19 * number ) >= AT.WorkTop ) {
733 MLOCK(ErrorMessageLock);
735 MUNLOCK(ErrorMessageLock);
742 for ( i = 0; i < number; i++ ) *p++ = *m++;
744 t->factor = params[4];
745 t->allsign = params[5];
746 if ( number >= 10 || ( t->gamma5 != GAMMA1 && number > 4 ) ) {
752 minimum = 0; min = t->inlist;
753 stopper = min + number;
754 for ( i = 1; i < number; i++ ) {
757 for ( j = 0; j < number; j++ ) {
758 if ( *p < *m )
break;
765 if ( m >= stopper ) m = t->inlist;
766 if ( p >= stopper ) p = t->inlist;
770 min = m = AT.WorkPointer;
774 if ( p >= stopper ) p = t->inlist;
779 while ( --i >= 0 ) *p++ = *m++;
783 i = *p; *p++ = *--m; *m = i;
786 for ( i = 0; i < number; i++ ) {
789 for ( j = 0; j < number; j++ ) {
790 if ( *p < *m )
break;
797 if ( m >= stopper ) m = t->inlist;
803 if ( m >= stopper ) m = t->inlist;
807 if ( ( minimum & 1 ) != 0 ) {
808 if ( t->gamma5 == GAMMA5 ) t->allsign = - t->allsign;
809 else if ( t->gamma5 != GAMMA1 )
810 t->gamma5 = GAMMA6 + GAMMA7 - t->gamma5;
812 p = min; m = t->inlist; i = number;
813 while ( --i >= 0 ) *m++ = *p++;
818 number = Trace4Gen(BHEAD t,number);
819 AT.WorkPointer = OldW;
856 WORD Trace4Gen(PHEAD
TRACES *t, WORD number)
859 WORD *termout, *stop;
861 WORD *pold, *mold, diff, *oldstring, cp;
866 if ( t->gamma5 == GAMMA5 )
return(0);
867 termout = AT.WorkPointer;
873 if ( *p == SUBEXPRESSION && p[2] == t->num ) {
876 do { *m++ = *p++; }
while ( p < oldstring );
878 *m++ = AC.lUniTrace[0];
879 *m++ = AC.lUniTrace[1];
880 *m++ = AC.lUniTrace[2];
881 *m++ = AC.lUniTrace[3];
882 if ( number == 2 || t->accup > t->accu ) {
890 if ( t->accup > t->accu ) {
893 while ( p < t->accup ) *m++ = *p++;
894 oldstring[1] = WORDDIF(m,oldstring);
904 do { *m++ = *p++; }
while ( p < stop );
905 *termout = WORDDIF(m,termout);
906 if ( t->allsign < 0 ) m[-1] = -m[-1];
907 if ( ( AT.WorkPointer = m ) > AT.WorkTop ) {
908 MLOCK(ErrorMessageLock);
910 MUNLOCK(ErrorMessageLock);
918 if (
Generator(BHEAD termout,t->level) )
goto TracCall;
920 AT.WorkPointer= termout;
924 }
while ( p < stop );
932 stop = p + number - 1;
938 while ( m < stop ) *p++ = *m++;
939 if ( t->gamma5 != GAMMA1 ) {
940 if ( t->gamma5 == GAMMA5 ) t->allsign = - t->allsign;
941 else if ( t->gamma5 == GAMMA6 ) t->gamma5 = GAMMA7;
942 else if ( t->gamma5 == GAMMA7 ) t->gamma5 = GAMMA6;
944 if ( Trace4Gen(BHEAD t,number-2) )
goto TracCall;
945 t = AN.tracestack + AN.numtracesctack - 1;
946 if ( t->gamma5 != GAMMA1 ) {
947 if ( t->gamma5 == GAMMA5 ) t->allsign = - t->allsign;
948 else if ( t->gamma5 == GAMMA6 ) t->gamma5 = GAMMA7;
949 else if ( t->gamma5 == GAMMA7 ) t->gamma5 = GAMMA6;
951 while ( p > t->inlist ) *--m = *--p;
963 while ( m <= stop ) *p++ = *m++;
964 if ( Trace4Gen(BHEAD t,number-2) )
goto TracCall;
965 t = AN.tracestack + AN.numtracesctack - 1;
966 while ( p > pold ) *--m = *--p;
973 }
while ( p < stop );
981 if ( *p >= AM.OffsetIndex && (
982 ( *p < WILDOFFSET + AM.OffsetIndex &&
983 indices[*p-AM.OffsetIndex].dimension == 4 )
984 || ( *p >= WILDOFFSET + AM.OffsetIndex && AC.lDefDim == 4 ) ) ) {
992 t->allsign = - t->allsign;
995 while ( m > p ) { diff = *p; *p++ = *m; *m-- = diff; }
998 while ( m < stop ) *p++ = *m++;
999 if ( Trace4Gen(BHEAD t,number-2) )
goto TracCall;
1000 t = AN.tracestack + AN.numtracesctack - 1;
1002 while ( m > mold ) *m-- = *--p;
1007 while ( m > p ) { diff = *p; *p++ = *m; *m-- = diff; }
1008 t->allsign = - t->allsign;
1016 }
while ( p < stop );
1025 if ( *p >= AM.OffsetIndex && (
1026 ( *p < WILDOFFSET + AM.OffsetIndex &&
1027 indices[*p-AM.OffsetIndex].dimension == 4 )
1028 || ( *p >= WILDOFFSET + AM.OffsetIndex && AC.lDefDim == 4 ) ) ) {
1030 if ( m >= stop ) m -= number;
1032 WORD oldfactor, old5;
1033 oldstring = AT.WorkPointer;
1034 AT.WorkPointer += number;
1035 oldfactor = t->allsign;
1037 if ( m < p ) cp = (WORDDIF(m,t->inlist) + 1 ) & 1;
1039 if ( cp && ( t->gamma5 != GAMMA1 ) ) {
1040 if ( t->gamma5 == GAMMA5 ) t->allsign = -t->allsign;
1041 else if ( t->gamma5 == GAMMA6 ) t->gamma5 = GAMMA7;
1042 else if ( t->gamma5 == GAMMA7 ) t->gamma5 = GAMMA6;
1047 while ( m < stop ) *p++ = *m++;
1053 while ( m < stop ) *p++ = *m++;
1055 if ( !cp && ((WORDDIF(stop,p))&1) != 0 && ( t->gamma5 != GAMMA1 ) ) {
1056 if ( t->gamma5 == GAMMA5 ) t->allsign = -t->allsign;
1057 else if ( t->gamma5 == GAMMA6 ) t->gamma5 = GAMMA7;
1058 else if ( t->gamma5 == GAMMA7 ) t->gamma5 = GAMMA6;
1060 while ( p < stop ) *p++ = *m++;
1064 oldval = number - 4;
1065 while ( oldval > 0 ) {
1066 if ( *p >= AM.OffsetIndex && (
1067 ( *p < WILDOFFSET + AM.OffsetIndex &&
1068 indices[*p-AM.OffsetIndex].dimension )
1069 || ( *p >= WILDOFFSET + AM.OffsetIndex && AC.lDefDim ) ) ) {
1074 else if ( *p == m[1] ) {
1081 if ( oldval <= 0 ) {
1082 *(t->accup)++ = *m++;
1083 *(t->accup)++ = *m++;
1085 if ( Trace4Gen(BHEAD t,number-4) )
goto TracCall;
1086 t = AN.tracestack + AN.numtracesctack - 1;
1088 if ( oldval <= 0 ) t->accup -= 2;
1090 t->allsign = oldfactor;
1091 AT.WorkPointer = p = oldstring;
1093 while ( m < stop ) *m++ = *p++;
1098 }
while ( p < stop );
1105 if ( *p >= AM.OffsetIndex && (
1106 ( *p < WILDOFFSET + AM.OffsetIndex &&
1107 indices[*p-AM.OffsetIndex].dimension == 4 )
1108 || ( *p >= WILDOFFSET + AM.OffsetIndex && AC.lDefDim == 4 ) ) ) {
1110 while ( m < stop ) {
1129 while ( m < stop ) *p++ = *m++;
1130 if ( Trace4Gen(BHEAD t,number-2) )
goto TracCall;
1131 t = AN.tracestack + AN.numtracesctack - 1;
1135 while ( m > p ) { diff = *p; *p++ = *m; *m-- = diff; }
1140 if ( Trace4Gen(BHEAD t,number-2) )
goto TracCall;
1141 t = AN.tracestack + AN.numtracesctack - 1;
1144 while ( m > p ) { diff = *p; *p++ = *m; *m-- = diff; }
1148 while ( m > mold ) *m-- = *--p;
1161 }
while ( p < stop );
1167 stop = p + number - 1;
1171 while ( p <= stop ) {
1173 if ( m > stop ) m -= number;
1175 WORD oldfactor, c, old5;
1176 oldfactor = t->allsign;
1178 cp = (WORDDIF(m,t->inlist)) & 1;
1179 if ( !cp && ( t->gamma5 != GAMMA1 ) ) {
1180 if ( t->gamma5 == GAMMA5 ) t->allsign = -t->allsign;
1181 else if ( t->gamma5 == GAMMA6 ) t->gamma5 = GAMMA7;
1182 else if ( t->gamma5 == GAMMA7 ) t->gamma5 = GAMMA6;
1184 oldstring = AT.WorkPointer;
1185 AT.WorkPointer += number;
1190 while ( m <= stop ) *p++ = *m++;
1196 while ( m <= stop ) *p++ = *m++;
1198 while ( p <= stop ) *p++ = *m++;
1202 *(t->accup) = oldval;
1209 if ( Trace4Gen(BHEAD t,number-2) )
goto Trac4Call;
1210 t = AN.tracestack + AN.numtracesctack - 1;
1212 t->allsign = - t->allsign;
1218 if ( Trace4Gen(BHEAD t,number-2) )
goto Trac4Call;
1219 t = AN.tracestack + AN.numtracesctack - 1;
1221 t->allsign = oldfactor;
1222 AT.WorkPointer = p = oldstring;
1224 while ( m <= stop ) *m++ = *p++;
1230 }
while ( ++diff <= (number>>1) );
1239 termout = AT.WorkPointer;
1241 if ( t->finalstep == 0 ) diff = Trace4no(number,t->accup,t);
1242 else diff = TraceNno(number,t->accup,t);
1244 if ( diff == 0 )
break;
1249 if ( p < stop )
do {
1250 if ( *p == SUBEXPRESSION && p[2] == t->num ) {
1253 do { *m++ = *p++; }
while ( p < oldstring );
1256 *m++ = AC.lUniTrace[0];
1257 *m++ = AC.lUniTrace[1];
1258 *m++ = AC.lUniTrace[2];
1259 *m++ = AC.lUniTrace[3];
1266 if ( diff == 2 || diff == -2 ) {
1271 *m++ = *p++; *m++ = *p++; *m++ = *p++; *m++ = *p++;
1274 if ( oldval > 0 || t->accup > t->accu ) {
1278 if ( oldval > 0 ) NCOPY(m,p,oldval);
1279 if ( t->accup > t->accu ) {
1281 while ( p < t->accup ) *m++ = *p++;
1282 oldstring[1] = WORDDIF(m,oldstring);
1286 do { *m++ = *p++; }
while ( p < stop );
1287 *termout = WORDDIF(m,termout);
1288 if ( ( diff ^ t->allsign ) < 0 ) m[-1] = - m[-1];
1289 if ( ( AT.WorkPointer = m ) > AT.WorkTop ) {
1290 MLOCK(ErrorMessageLock);
1292 MUNLOCK(ErrorMessageLock);
1298 if (
Generator(BHEAD termout,t->level) ) {
1299 AT.WorkPointer = termout;
1302 t = AN.tracestack + AN.numtracesctack - 1;
1307 }
while ( p < stop );
1309 AT.WorkPointer = termout;
1316 AT.WorkPointer = oldstring;
1318 if ( AM.tracebackflag ) {
1319 MLOCK(ErrorMessageLock);
1320 MesCall(
"Trace4Gen");
1321 MUNLOCK(ErrorMessageLock);
1341 WORD TraceNno(WORD number, WORD *kron,
TRACES *t)
1345 if ( !number || ( number & 1 ) )
return(0);
1347 for ( i = 0; i < number; i++ ) {
1359 for ( j = i + 1; j <= *p; j++ ) kron[j-1] = kron[j];
1361 if ( *p < number ) {
1364 while ( j >= (i+1) ) { kron[j] = kron[j-1]; j--; }
1367 for ( j = i+2; j < number; j += 2 ) t->perm[j] = j;
1382 WORD TraceN(PHEAD WORD *term, WORD *params, WORD num, WORD level)
1386 WORD *p, *m, number, i;
1388 if ( params[3] != GAMMA1 ) {
1389 MLOCK(ErrorMessageLock);
1390 MesPrint(
"Gamma5 not allowed in n-trace");
1391 MUNLOCK(ErrorMessageLock);
1394 OldW = AT.WorkPointer;
1395 if ( AN.numtracesctack >= AN.intracestack ) {
1396 number = AN.intracestack + 2;
1397 t = (
TRACES *)Malloc1(number*
sizeof(
TRACES),
"TRACES-struct");
1398 if ( AN.tracestack ) {
1399 for ( i = 0; i < AN.intracestack; i++ ) { t[i] = AN.tracestack[i]; }
1400 M_free(AN.tracestack,
"TRACES-struct");
1403 AN.intracestack = number;
1405 number = *params - 6;
1406 if ( number < 0 || ( number & 1 ) || !params[5] )
return(0);
1408 t = AN.tracestack + AN.numtracesctack;
1409 AN.numtracesctack++;
1411 t->inlist = AT.WorkPointer;
1412 t->accup = t->accu = t->inlist + number;
1413 t->perm = t->accu + number;
1414 if ( ( AT.WorkPointer += 3 * number ) >= AT.WorkTop ) {
1415 AN.numtracesctack--;
1416 MLOCK(ErrorMessageLock);
1418 MUNLOCK(ErrorMessageLock);
1425 for ( i = 0; i < number; i++ ) *p++ = *m++;
1427 t->factor = params[4];
1428 t->allsign = params[5];
1429 number = TraceNgen(BHEAD t,number);
1430 AT.WorkPointer = OldW;
1431 AN.numtracesctack--;
1446 WORD TraceNgen(PHEAD
TRACES *t, WORD number)
1449 WORD *termout, *stop;
1450 WORD *p, *m, oldval;
1451 WORD *pold, *mold, diff, *oldstring;
1455 if ( number <= 2 ) {
1456 termout = AT.WorkPointer;
1461 if ( p < stop )
do {
1462 if ( *p == SUBEXPRESSION && p[2] == t->num ) {
1465 do { *m++ = *p++; }
while ( p < oldstring );
1467 *m++ = AC.lUniTrace[0];
1468 *m++ = AC.lUniTrace[1];
1469 *m++ = AC.lUniTrace[2];
1470 *m++ = AC.lUniTrace[3];
1471 if ( number == 2 || t->accup > t->accu ) {
1475 if ( number == 2 ) {
1476 *m++ = t->inlist[0];
1477 *m++ = t->inlist[1];
1479 if ( t->accup > t->accu ) {
1482 while ( p < t->accup ) *m++ = *p++;
1483 oldstring[1] = WORDDIF(m,oldstring);
1493 do { *m++ = *p++; }
while ( p < stop );
1494 *termout = WORDDIF(m,termout);
1495 if ( t->allsign < 0 ) m[-1] = -m[-1];
1496 if ( ( AT.WorkPointer = m ) > AT.WorkTop ) {
1497 MLOCK(ErrorMessageLock);
1499 MUNLOCK(ErrorMessageLock);
1505 if (
Generator(BHEAD termout,t->level) )
goto TracCall;
1507 AT.WorkPointer= termout;
1511 }
while ( p < stop );
1519 stop = p + number - 1;
1520 if ( *p == *stop ) {
1525 while ( m < stop ) *p++ = *m++;
1526 if ( TraceNgen(BHEAD t,number-2) )
goto TracCall;
1527 t = AN.tracestack + AN.numtracesctack - 1;
1528 while ( p > t->inlist ) *--m = *--p;
1529 *p = *stop = oldval;
1540 while ( m <= stop ) *p++ = *m++;
1541 if ( TraceNgen(BHEAD t,number-2) )
goto TracCall;
1542 t = AN.tracestack + AN.numtracesctack - 1;
1543 while ( p > pold ) *--m = *--p;
1550 }
while ( p < stop );
1556 stop = p + number - 1;
1560 while ( p <= stop ) {
1562 if ( m > stop ) m -= number;
1565 oldstring = AT.WorkPointer;
1566 AT.WorkPointer += number;
1571 while ( m <= stop ) *p++ = *m++;
1578 while ( m <= stop ) *p++ = *m++;
1580 while ( p <= stop ) *p++ = *m++;
1582 oldfactor = t->allsign;
1586 if ( oldval >= ( AM.OffsetIndex + WILDOFFSET ) ||
1587 ( oldval >= AM.OffsetIndex
1588 && indices[oldval-AM.OffsetIndex].dimension ) ) {
1598 while ( m > (p+3) ) {
1602 if ( TraceNgen(BHEAD t,number-2) )
goto TracnCall;
1603 t = AN.tracestack + AN.numtracesctack - 1;
1605 t->allsign = - t->allsign;
1607 switch ( WORDDIF(m,p) ) {
1612 if ( oldval < ( AM.OffsetIndex + WILDOFFSET )
1613 && indices[oldval-AM.OffsetIndex].nmin4
1615 t->allsign = - t->allsign;
1616 if ( TraceNgen(BHEAD t,number-2) )
goto TracnCall;
1617 t = AN.tracestack + AN.numtracesctack - 1;
1619 *(t->accup)++ = SUMMEDIND;
1621 indices[oldval-AM.OffsetIndex].nmin4;
1625 if ( TraceNgen(BHEAD t,number-2) )
goto TracnCall;
1626 t = AN.tracestack + AN.numtracesctack - 1;
1627 t->allsign = - t->allsign;
1629 *(t->accup)++ = oldval;
1630 *(t->accup)++ = oldval;
1632 if ( TraceNgen(BHEAD t,number-2) )
goto TracnCall;
1633 t = AN.tracestack + AN.numtracesctack - 1;
1643 if ( TraceNgen(BHEAD t,number-4) )
goto TracnCall;
1644 t = AN.tracestack + AN.numtracesctack - 1;
1645 *p = one; p[1] = two;
1647 if ( oldval < ( AM.OffsetIndex + WILDOFFSET )
1648 && indices[oldval-AM.OffsetIndex].nmin4
1651 *(t->accup)++ = SUMMEDIND;
1653 indices[oldval-AM.OffsetIndex].nmin4;
1656 t->allsign = - t->allsign;
1657 if ( TraceNgen(BHEAD t,number-2) )
goto TracnCall;
1658 t = AN.tracestack + AN.numtracesctack - 1;
1659 t->allsign = - t->allsign;
1661 *(t->accup)++ = oldval;
1662 *(t->accup)++ = oldval;
1664 if ( TraceNgen(BHEAD t,number-2) )
goto TracnCall;
1665 t = AN.tracestack + AN.numtracesctack - 1;
1673 c = m[-1]; m[-1] = m[-2]; m[-2] = c;
1674 t->allsign = - t->allsign;
1675 if ( TraceNgen(BHEAD t,number-2) )
goto TracnCall;
1676 t = AN.tracestack + AN.numtracesctack - 1;
1682 if ( oldval < ( AM.OffsetIndex + WILDOFFSET )
1683 && indices[oldval-AM.OffsetIndex].nmin4
1685 *(t->accup)++ = SUMMEDIND;
1687 indices[oldval-AM.OffsetIndex].nmin4;
1688 if ( TraceNgen(BHEAD t,number-2) )
goto TracnCall;
1689 t = AN.tracestack + AN.numtracesctack - 1;
1691 t->allsign = - t->allsign;
1695 *(t->accup)++ = oldval;
1696 *(t->accup)++ = oldval;
1697 if ( TraceNgen(BHEAD t,number-2) )
goto TracnCall;
1698 t = AN.tracestack + AN.numtracesctack - 1;
1700 t->allsign = - t->allsign;
1702 if ( TraceNgen(BHEAD t,number-2) )
goto TracnCall;
1703 t = AN.tracestack + AN.numtracesctack - 1;
1710 *(t->accup) = oldval;
1717 if ( TraceNgen(BHEAD t,number-2) )
goto TracnCall;
1718 t = AN.tracestack + AN.numtracesctack - 1;
1720 t->allsign = - t->allsign;
1726 if ( TraceNgen(BHEAD t,number-2) )
goto TracnCall;
1727 t = AN.tracestack + AN.numtracesctack - 1;
1730 t->allsign = oldfactor;
1733 while ( m <= stop ) *m++ = *p++;
1734 AT.WorkPointer = oldstring;
1740 }
while ( diff <= (number>>1) );
1749 termout = AT.WorkPointer;
1750 while ( ( diff = TraceNno(number,t->accup,t) ) != 0 ) {
1755 if ( p < stop )
do {
1756 if ( *p == SUBEXPRESSION && p[2] == t->num ) {
1759 do { *m++ = *p++; }
while ( p < oldstring );
1762 *m++ = AC.lUniTrace[0];
1763 *m++ = AC.lUniTrace[1];
1764 *m++ = AC.lUniTrace[2];
1765 *m++ = AC.lUniTrace[3];
1776 if ( t->accup > t->accu ) {
1778 while ( p < t->accup ) *m++ = *p++;
1779 oldstring[1] = WORDDIF(m,oldstring);
1782 do { *m++ = *p++; }
while ( p < stop );
1783 *termout = WORDDIF(m,termout);
1784 if ( ( diff ^ t->allsign ) < 0 ) m[-1] = - m[-1];
1785 if ( ( AT.WorkPointer = m ) > AT.WorkTop ) {
1786 MLOCK(ErrorMessageLock);
1788 MUNLOCK(ErrorMessageLock);
1794 if (
Generator(BHEAD termout,t->level) ) {
1795 AT.WorkPointer = termout;
1798 t = AN.tracestack + AN.numtracesctack - 1;
1803 }
while ( p < stop );
1805 AT.WorkPointer = termout;
1812 AT.WorkPointer = oldstring;
1814 if ( AM.tracebackflag ) {
1815 MLOCK(ErrorMessageLock);
1816 MesCall(
"TraceNGen");
1817 MUNLOCK(ErrorMessageLock);
1831 WORD Traces(PHEAD WORD *term, WORD *params, WORD num, WORD level)
1834 switch ( AT.TMout[2] ) {
1836 return(TraceN(BHEAD term,params,num,level));
1838 return(Trace4(BHEAD term,params,num,level));
1840 return(Trace4(BHEAD term,params,num,level));
1842 return(Trace4(BHEAD term,params,num,level));
1853 WORD TraceFind(PHEAD WORD *term, WORD *params)
1857 WORD *termout, *stop, *stop2, number = 0;
1859 WORD type, spinline, sp;
1861 spinline = params[4];
1862 if ( spinline < 0 ) {
1863 sp = DolToIndex(BHEAD -spinline);
1864 if ( AN.ErrorInDollar || sp < 0 ) {
1865 MLOCK(ErrorMessageLock);
1866 MesPrint(
"$%s does not have an index value in trace statement in module %l",
1867 DOLLARNAME(Dollars,-spinline),AC.CModule);
1868 MUNLOCK(ErrorMessageLock);
1883 termout = m = AT.WorkPointer;
1886 while ( p < stop ) {
1888 if ( *p == GAMMA && p[FUNHEAD] == spinline ) {
1890 *m++ = SUBEXPRESSION;
1899 while ( p < stop2 ) {
1900 if ( *p == GAMMA5 ) {
1901 if ( AT.TMout[3] == GAMMA5 ) AT.TMout[3] = GAMMA1;
1902 else if ( AT.TMout[3] == GAMMA1 ) AT.TMout[3] = GAMMA5;
1903 else if ( AT.TMout[3] == GAMMA7 ) AT.TMout[5] = -AT.TMout[5];
1904 if ( number & 1 ) AT.TMout[5] = - AT.TMout[5];
1907 else if ( *p == GAMMA6 ) {
1908 if ( number & 1 )
goto F7;
1909 F6:
if ( AT.TMout[3] == GAMMA6 ) (AT.TMout[4])++;
1910 else if ( AT.TMout[3] == GAMMA1 ) AT.TMout[3] = GAMMA6;
1911 else if ( AT.TMout[3] == GAMMA5 ) AT.TMout[3] = GAMMA6;
1912 else if ( AT.TMout[3] == GAMMA7 ) AT.TMout[5] = 0;
1915 else if ( *p == GAMMA7 ) {
1916 if ( number & 1 )
goto F6;
1917 F7:
if ( AT.TMout[3] == GAMMA7 ) (AT.TMout[4])++;
1918 else if ( AT.TMout[3] == GAMMA1 ) AT.TMout[3] = GAMMA7;
1919 else if ( AT.TMout[3] == GAMMA5 ) {
1920 AT.TMout[3] = GAMMA7;
1921 AT.TMout[5] = -AT.TMout[5];
1923 else if ( AT.TMout[3] == GAMMA6 ) AT.TMout[5] = 0;
1933 while ( p < stop2 ) *m++ = *p++;
1936 if ( first )
return(0);
1937 AT.TMout[0] = WORDDIF(to,AT.TMout);
1940 while ( p < to ) *m++ = *p++;
1941 *termout = WORDDIF(m,termout);
1944 do { *to++ = *p++; }
while ( p < m );
1945 AT.WorkPointer = term + *term;
1963 WORD Chisholm(PHEAD WORD *term, WORD level)
1966 WORD *t, *r, *m, *s, *tt, *rr;
1967 WORD *mat, *matpoint, *termout, *rdo;
1968 CBUF *C = cbuf+AM.rbufnum;
1969 WORD i, j, num = C->
lhs[level][2], gam5;
1970 WORD norm = 0, k, *matp;
1974 mat = matpoint = AT.WorkPointer;
1976 r = t + *t - 1; r -= ABS(*r);
1981 if ( *t == GAMMA && t[FUNHEAD] == num ) {
1985 if ( *t >= 0 || *t < MINSPEC ) i++;
1987 if ( gam5 == GAMMA1 ) gam5 = *t;
1988 else if ( gam5 == GAMMA5 ) {
1989 if ( *t == GAMMA5 ) gam5 = GAMMA1;
1990 else if ( *t != GAMMA1 ) gam5 = *t;
1998 if ( ( i & 1 ) != 0 )
return(0);
2014 while ( s < matpoint ) {
2019 if ( *s < AM.OffsetIndex || ( *s < ( AM.OffsetIndex + WILDOFFSET ) &&
2020 indices[*s-AM.OffsetIndex].dimension != 4 )
2021 || ( ( AC.lDefDim != 4 ) && ( *s >= ( AM.OffsetIndex + WILDOFFSET ) ) ) ) {
2026 if ( *t == GAMMA && t[FUNHEAD] != num ) {
2040 if ( norm == 0 )
return(
Generator(BHEAD term,level));
2060 if ( C->
lhs[level][3] == 0 ) norm = 1;
2063 for ( k = 0; k < norm; k++ ) {
2066 while ( s < matpoint ) {
2071 if ( *s < AM.OffsetIndex || ( *s < ( AM.OffsetIndex + WILDOFFSET ) &&
2072 indices[*s-AM.OffsetIndex].dimension != 4 ) ) {
2077 if ( *t == GAMMA && t[FUNHEAD] != num ) {
2085 while ( m <= s ) *matpoint++ = *m++;
2087 while ( m < matpoint ) *t++ = *m++;
2092 if ( *t != GAMMA || t[FUNHEAD] != num ) {
2102 while ( --j >= 0 ) *m++ = *t++;
2105 while ( s < termout ) *m++ = *s++;
2108 while ( t < tt ) *m++ = *t++;
2109 rdo[1] = WORDDIF(m,rdo);
2111 *m++ = AC.lUniTrace[0];
2112 *m++ = AC.lUniTrace[1];
2113 *m++ = AC.lUniTrace[2];
2114 *m++ = AC.lUniTrace[3];
2121 if ( *t != GAMMA || t[FUNHEAD] != num ) {
2128 while ( t < rr ) *m++ = *t++;
2130 *termout = WORDDIF(m,termout);
2136 if (
Generator(BHEAD t,level) )
goto ChisCall;
2138 j = WORDDIF(termout,mat)-1;
2141 AT.WorkPointer = rr;
2143 i = *--m; *m = *t; *t++ = i;
2146 if (
Generator(BHEAD termout,level) )
goto ChisCall;
2147 AT.WorkPointer = mat;
2165 if ( AM.tracebackflag ) {
2166 MLOCK(ErrorMessageLock);
2167 MesCall(
"Chisholm");
2168 MUNLOCK(ErrorMessageLock);
2178 WORD TenVecFind(PHEAD WORD *term, WORD *params)
2181 WORD *t, *w, *m, *tstop;
2182 WORD i, mode, thevector, thetensor, spectator;
2183 thetensor = params[3];
2184 thevector = params[4];
2186 if ( thetensor < 0 ) {
2187 thetensor = DolToTensor(BHEAD -thetensor);
2188 if ( thetensor < FUNCTION ) {
2189 if ( thevector > 0 ) {
2190 thetensor = DolToTensor(BHEAD thevector);
2191 if ( thetensor < FUNCTION ) {
2192 MLOCK(ErrorMessageLock);
2193 MesPrint(
"$%s should have been a tensor in module %l"
2194 ,DOLLARNAME(Dollars,params[4]),AC.CModule);
2195 MUNLOCK(ErrorMessageLock);
2198 thevector = DolToVector(BHEAD -params[3]);
2199 if ( thevector >= 0 ) {
2200 MLOCK(ErrorMessageLock);
2201 MesPrint(
"$%s should have been a vector in module %l"
2202 ,DOLLARNAME(Dollars,-params[3]),AC.CModule);
2203 MUNLOCK(ErrorMessageLock);
2208 MLOCK(ErrorMessageLock);
2209 MesPrint(
"$%s should have been a tensor in module %l"
2210 ,DOLLARNAME(Dollars,-params[3]),AC.CModule);
2211 MUNLOCK(ErrorMessageLock);
2216 if ( thevector > 0 ) {
2217 thevector = DolToVector(BHEAD thevector);
2218 if ( thevector >= 0 ) {
2219 MLOCK(ErrorMessageLock);
2220 MesPrint(
"$%s should have been a vector in module %l"
2221 ,DOLLARNAME(Dollars,params[4]),AC.CModule);
2222 MUNLOCK(ErrorMessageLock);
2226 if ( ( mode & 1 ) != 0 ) {
2227 GETSTOP(term,tstop);
2229 while ( t < tstop ) {
2230 if ( *t == DOTPRODUCT ) {
2231 i = t[1] - 2; t += 2;
2235 else if ( *t == thevector && t[1] == thevector ) {
2236 if ( ( mode & 2 ) == 0 ) spectator = thevector;
2238 else if ( *t == thevector ) spectator = t[1];
2239 else if ( t[1] == thevector ) spectator = *t;
2241 if ( ( mode & 8 ) == 0 )
goto match;
2242 w = SetElements + Sets[params[6]].first;
2243 m = SetElements + Sets[params[6]].last;
2245 if ( *w == spectator )
break;
2248 if ( w >= m )
goto match;
2254 else if ( *t == VECTOR ) {
2255 i = t[1] - 2; t += 2;
2257 if ( *t == thevector )
goto match;
2262 else if ( *t == thetensor ) t += t[1];
2263 else if ( *t >= FUNCTION ) {
2264 if ( functions[*t-FUNCTION].spec > 0 ) {
2268 if ( *t == thevector )
goto match;
2272 else if ( ( mode & 4 ) != 0 ) {
2276 if ( *t == -VECTOR && t[1] == thevector )
goto match;
2277 else if ( *t > 0 ) t += *t;
2278 else if ( *t <= -FUNCTION ) t++;
2288 GETSTOP(term,tstop);
2290 while ( t < tstop ) {
2291 if ( *t == thetensor )
goto match;
2298 AT.TMout[1] = TENVEC;
2299 AT.TMout[2] = thetensor;
2300 AT.TMout[3] = thevector;
2302 if ( ( mode & 8 ) != 0 ) { AT.TMout[0] = 6; AT.TMout[5] = params[6]; }
2312 WORD TenVec(PHEAD WORD *term, WORD *params, WORD num, WORD level)
2315 WORD *t, *m, *w, *termout, *tstop, *outlist, *ou, *ww, *mm;
2316 WORD i, j, k, x, mode, thevector, thetensor, DumNow, spectator;
2318 thetensor = params[2];
2319 thevector = params[3];
2321 termout = AT.WorkPointer;
2322 DumNow = AR.CurDum = DetCurDum(BHEAD term);
2323 if ( ( mode & 1 ) != 0 ) {
2324 AT.WorkPointer += *term;
2325 ou = outlist = AT.WorkPointer;
2326 GETSTOP(term,tstop);
2329 while ( t < tstop ) {
2330 if ( *t == DOTPRODUCT ) {
2333 *m++ = *t++; *m++ = *t++;
2337 *m++ = *t++; *m++ = *t++; *m++ = *t++;
2339 else if ( *t == thevector && t[1] == thevector ) {
2340 if ( ( mode & 2 ) == 0 ) spectator = thevector;
2342 *m++ = *t++; *m++ = *t++; *m++ = *t++;
2345 else if ( *t == thevector ) spectator = t[1];
2346 else if ( t[1] == thevector ) spectator = *t;
2348 *m++ = *t++; *m++ = *t++; *m++ = *t++;
2351 if ( ( mode & 8 ) == 0 )
goto noveto;
2352 ww = SetElements + Sets[params[5]].first;
2353 mm = SetElements + Sets[params[5]].last;
2355 if ( *ww == spectator )
break;
2359 *m++ = *t++; *m++ = *t++; *m++ = *t++;
2362 noveto:
if ( spectator == thevector ) {
2363 for ( j = 0; j < t[2]; j++ ) {
2364 *ou++ = ++AR.CurDum;
2370 for ( j = 0; j < t[2]; j++ ) *ou++ = spectator;
2376 w[1] = WORDDIF(m,w);
2377 if ( w[1] == 2 ) m = w;
2379 else if ( *t == VECTOR ) {
2380 i = t[1] - 2; w = m;
2381 *m++ = *t++; *m++ = *t++;
2383 if ( *t == thevector ) {
2387 else { *m++ = *t++; *m++ = *t++; }
2390 w[1] = WORDDIF(m,w);
2391 if ( w[1] == 2 ) m = w;
2393 else if ( *t == thetensor ) {
2398 else if ( *t >= FUNCTION ) {
2399 if ( functions[*t-FUNCTION].spec > 0 ) {
2404 if ( *t == thevector ) {
2412 else if ( ( mode & 4 ) != 0 ) {
2417 if ( *t == -VECTOR && t[1] == thevector ) {
2423 else if ( *t > 0 ) {
2427 else if ( *t <= -FUNCTION ) *m++ = *t++;
2428 else { *m++ = *t++; *m++ = *t++; }
2439 i = WORDDIF(ou,outlist);
2441 for ( j = 1; j < i; j++ ) {
2442 if ( outlist[j-1] > outlist[j] ) {
2443 x = outlist[j-1]; outlist[j-1] = outlist[j]; outlist[j] = x;
2444 for ( k = j-1; k > 0; k-- ) {
2445 if ( outlist[k-1] <= outlist[k] )
break;
2446 x = outlist[k-1]; outlist[k-1] = outlist[k]; outlist[k] = x;
2453 *m++ = DIRTYSYMFLAG;
2459 while ( t < w ) *m++ = *t++;
2462 GETSTOP(term,tstop);
2465 while ( t < tstop ) {
2466 if ( *t != thetensor ) {
2475 while ( --i >= 0 ) {
2480 w[1] = WORDDIF(m,w);
2485 while ( t < w ) *m++ = *t++;
2487 *termout = WORDDIF(m,termout);
2490 if (
Generator(BHEAD termout,level) )
goto fromTenVec;
2492 AT.WorkPointer = termout;
2495 if ( AM.tracebackflag ) {
2496 MLOCK(ErrorMessageLock);
2498 MUNLOCK(ErrorMessageLock);
WORD Generator(PHEAD WORD *, WORD)