[6d07e3] | 1 | /***************************************** |
---|
| 2 | * Computer Algebra System SINGULAR * |
---|
| 3 | *****************************************/ |
---|
| 4 | /* |
---|
| 5 | * ABSTRACT: Gauss-Manin system normal form |
---|
| 6 | */ |
---|
| 7 | |
---|
[9f7665] | 8 | |
---|
| 9 | |
---|
| 10 | |
---|
[b1dfaf] | 11 | #include <kernel/mod2.h> |
---|
[6d07e3] | 12 | |
---|
| 13 | #ifdef HAVE_GMS |
---|
| 14 | |
---|
[4dbeb1] | 15 | #include "gms.h" |
---|
| 16 | |
---|
[0fb34ba] | 17 | #include <coeffs/numbers.h> |
---|
[737a68] | 18 | #include <kernel/polys.h> |
---|
[6d07e3] | 19 | |
---|
[4dbeb1] | 20 | #include "ipid.h" |
---|
[a7c2db] | 21 | |
---|
[5ca9807] | 22 | lists gmsNF(ideal p,ideal g,matrix B,int D,int K) |
---|
[6d07e3] | 23 | { |
---|
| 24 | ideal r=idInit(IDELEMS(p),1); |
---|
| 25 | ideal q=idInit(IDELEMS(p),1); |
---|
| 26 | |
---|
[a7c2db] | 27 | matrix B0=mpNew(MATROWS(B),MATCOLS(B)); |
---|
| 28 | for(int i=1;i<=MATROWS(B0);i++) |
---|
[4dbeb1] | 29 | for(int j=1;j<=MATCOLS(B0);j++) |
---|
[a7c2db] | 30 | if(MATELEM(B,i,j)!=NULL) |
---|
| 31 | MATELEM(B0,i,j)=pDiff(MATELEM(B,i,j),i+1); |
---|
[6d07e3] | 32 | |
---|
[a7c2db] | 33 | for(int k=0;k<IDELEMS(p);k++) |
---|
[6d07e3] | 34 | { |
---|
| 35 | while(p->m[k]!=NULL&&pGetExp(p->m[k],1)<=K) |
---|
| 36 | { |
---|
[a7c2db] | 37 | int j=0; |
---|
[6d07e3] | 38 | while(j<IDELEMS(g)&&!pLmDivisibleBy(g->m[j],p->m[k])) |
---|
| 39 | j++; |
---|
| 40 | |
---|
| 41 | if(j<IDELEMS(g)) |
---|
| 42 | { |
---|
[a7c2db] | 43 | poly m=pDivideM(pHead(p->m[k]),pHead(g->m[j])); |
---|
| 44 | p->m[k]=pSub(p->m[k],ppMult_mm(g->m[j],m)); |
---|
| 45 | pIncrExp(m,1); |
---|
| 46 | pSetm(m); |
---|
| 47 | for(int i=0;i<MATROWS(B);i++) |
---|
[18f7ec] | 48 | { |
---|
[a7c2db] | 49 | poly m0=pDiff(m,i+2); |
---|
| 50 | if(MATELEM(B0,i+1,j+1)!=NULL) |
---|
| 51 | p->m[k]=pAdd(p->m[k],ppMult_mm(MATELEM(B0,i+1,j+1),m)); |
---|
| 52 | if(MATELEM(B,i+1,j+1)!=NULL&&m0!=NULL) |
---|
| 53 | p->m[k]=pAdd(p->m[k],ppMult_mm(MATELEM(B,i+1,j+1),m0)); |
---|
| 54 | pDelete(&m0); |
---|
| 55 | } |
---|
| 56 | pDelete(&m); |
---|
[6d07e3] | 57 | } |
---|
| 58 | else |
---|
| 59 | { |
---|
[a7c2db] | 60 | poly p0=p->m[k]; |
---|
| 61 | pIter(p->m[k]); |
---|
| 62 | pNext(p0)=NULL; |
---|
| 63 | r->m[k]=pAdd(r->m[k],p0); |
---|
[6d07e3] | 64 | } |
---|
| 65 | |
---|
[a7c2db] | 66 | while(p->m[k]!=NULL&&pGetExp(p->m[k],1)<=K&&pWTotaldegree(p->m[k])>D) |
---|
[6d07e3] | 67 | { |
---|
[a7c2db] | 68 | int i=pGetExp(p->m[k],1); |
---|
| 69 | do |
---|
[18f7ec] | 70 | { |
---|
[a7c2db] | 71 | poly p0=p->m[k]; |
---|
| 72 | pIter(p->m[k]); |
---|
| 73 | pNext(p0)=NULL; |
---|
| 74 | q->m[k]=pAdd(q->m[k],p0); |
---|
| 75 | }while(p->m[k]!=NULL&&pGetExp(p->m[k],1)==i); |
---|
[6d07e3] | 76 | } |
---|
[a7c2db] | 77 | |
---|
| 78 | pNormalize(p->m[k]); |
---|
[6d07e3] | 79 | } |
---|
| 80 | |
---|
[a7c2db] | 81 | q->m[k]=pAdd(q->m[k],p->m[k]); |
---|
| 82 | p->m[k]=NULL; |
---|
[6d07e3] | 83 | } |
---|
[a7c2db] | 84 | idDelete(&p); |
---|
| 85 | idDelete((ideal *)&B0); |
---|
| 86 | |
---|
[4dbeb1] | 87 | id_Normalize(r, currRing); |
---|
| 88 | id_Normalize(q, currRing); |
---|
[a7c2db] | 89 | |
---|
| 90 | lists l=(lists)omAllocBin(slists_bin); |
---|
| 91 | l->Init(2); |
---|
| 92 | |
---|
| 93 | l->m[0].rtyp=IDEAL_CMD; |
---|
| 94 | l->m[0].data=r; |
---|
| 95 | l->m[1].rtyp=IDEAL_CMD; |
---|
| 96 | l->m[1].data=q; |
---|
[6d07e3] | 97 | |
---|
| 98 | return l; |
---|
| 99 | } |
---|
| 100 | |
---|
[a7c2db] | 101 | |
---|
[5ca9807] | 102 | BOOLEAN gmsNF(leftv res,leftv h) |
---|
[6d07e3] | 103 | { |
---|
| 104 | if(currRingHdl) |
---|
| 105 | { |
---|
| 106 | if(h&&h->Typ()==IDEAL_CMD) |
---|
| 107 | { |
---|
[18f7ec] | 108 | ideal p=(ideal)h->CopyD(); |
---|
[6d07e3] | 109 | h=h->next; |
---|
| 110 | if(h&&h->Typ()==IDEAL_CMD) |
---|
| 111 | { |
---|
| 112 | ideal g=(ideal)h->Data(); |
---|
| 113 | h=h->next; |
---|
| 114 | if(h&&h->Typ()==MATRIX_CMD) |
---|
| 115 | { |
---|
| 116 | matrix B=(matrix)h->Data(); |
---|
| 117 | h=h->next; |
---|
| 118 | if(h&&h->Typ()==INT_CMD) |
---|
| 119 | { |
---|
[6caadfb] | 120 | int D=(int)(long)h->Data(); |
---|
[6d07e3] | 121 | h=h->next; |
---|
| 122 | if(h&&h->Typ()==INT_CMD) |
---|
| 123 | { |
---|
[6caadfb] | 124 | int K=(int)(long)h->Data(); |
---|
[6d07e3] | 125 | res->rtyp=LIST_CMD; |
---|
[18f7ec] | 126 | res->data=(void *)gmsNF(p,g,B,D,K); |
---|
[6d07e3] | 127 | return FALSE; |
---|
[18f7ec] | 128 | } |
---|
| 129 | } |
---|
[6d07e3] | 130 | } |
---|
| 131 | } |
---|
| 132 | } |
---|
| 133 | WerrorS("<ideal>,<ideal>,<matrix>,<int>,<int> expected"); |
---|
| 134 | return TRUE; |
---|
| 135 | } |
---|
| 136 | WerrorS("no ring active"); |
---|
| 137 | return TRUE; |
---|
| 138 | } |
---|
| 139 | |
---|
| 140 | #endif /* HAVE_GMS */ |
---|
[398ad5] | 141 | |
---|