[BACK]Return to mat.c CVS log [TXT][DIR] Up to [local] / OpenXM_contrib2 / asir2000 / engine

Annotation of OpenXM_contrib2/asir2000/engine/mat.c, Revision 1.2

1.2     ! noro        1: /*
        !             2:  * Copyright (c) 1994-2000 FUJITSU LABORATORIES LIMITED
        !             3:  * All rights reserved.
        !             4:  *
        !             5:  * FUJITSU LABORATORIES LIMITED ("FLL") hereby grants you a limited,
        !             6:  * non-exclusive and royalty-free license to use, copy, modify and
        !             7:  * redistribute, solely for non-commercial and non-profit purposes, the
        !             8:  * computer program, "Risa/Asir" ("SOFTWARE"), subject to the terms and
        !             9:  * conditions of this Agreement. For the avoidance of doubt, you acquire
        !            10:  * only a limited right to use the SOFTWARE hereunder, and FLL or any
        !            11:  * third party developer retains all rights, including but not limited to
        !            12:  * copyrights, in and to the SOFTWARE.
        !            13:  *
        !            14:  * (1) FLL does not grant you a license in any way for commercial
        !            15:  * purposes. You may use the SOFTWARE only for non-commercial and
        !            16:  * non-profit purposes only, such as academic, research and internal
        !            17:  * business use.
        !            18:  * (2) The SOFTWARE is protected by the Copyright Law of Japan and
        !            19:  * international copyright treaties. If you make copies of the SOFTWARE,
        !            20:  * with or without modification, as permitted hereunder, you shall affix
        !            21:  * to all such copies of the SOFTWARE the above copyright notice.
        !            22:  * (3) An explicit reference to this SOFTWARE and its copyright owner
        !            23:  * shall be made on your publication or presentation in any form of the
        !            24:  * results obtained by use of the SOFTWARE.
        !            25:  * (4) In the event that you modify the SOFTWARE, you shall notify FLL by
        !            26:  * e-mail at risa-admin@flab.fujitsu.co.jp of the detailed specification
        !            27:  * for such modification or the source code of the modified part of the
        !            28:  * SOFTWARE.
        !            29:  *
        !            30:  * THE SOFTWARE IS PROVIDED AS IS WITHOUT ANY WARRANTY OF ANY KIND. FLL
        !            31:  * MAKES ABSOLUTELY NO WARRANTIES, EXPRESSED, IMPLIED OR STATUTORY, AND
        !            32:  * EXPRESSLY DISCLAIMS ANY IMPLIED WARRANTY OF MERCHANTABILITY, FITNESS
        !            33:  * FOR A PARTICULAR PURPOSE OR NONINFRINGEMENT OF THIRD PARTIES'
        !            34:  * RIGHTS. NO FLL DEALER, AGENT, EMPLOYEES IS AUTHORIZED TO MAKE ANY
        !            35:  * MODIFICATIONS, EXTENSIONS, OR ADDITIONS TO THIS WARRANTY.
        !            36:  * UNDER NO CIRCUMSTANCES AND UNDER NO LEGAL THEORY, TORT, CONTRACT,
        !            37:  * OR OTHERWISE, SHALL FLL BE LIABLE TO YOU OR ANY OTHER PERSON FOR ANY
        !            38:  * DIRECT, INDIRECT, SPECIAL, INCIDENTAL, PUNITIVE OR CONSEQUENTIAL
        !            39:  * DAMAGES OF ANY CHARACTER, INCLUDING, WITHOUT LIMITATION, DAMAGES
        !            40:  * ARISING OUT OF OR RELATING TO THE SOFTWARE OR THIS AGREEMENT, DAMAGES
        !            41:  * FOR LOSS OF GOODWILL, WORK STOPPAGE, OR LOSS OF DATA, OR FOR ANY
        !            42:  * DAMAGES, EVEN IF FLL SHALL HAVE BEEN INFORMED OF THE POSSIBILITY OF
        !            43:  * SUCH DAMAGES, OR FOR ANY CLAIM BY ANY OTHER PARTY. EVEN IF A PART
        !            44:  * OF THE SOFTWARE HAS BEEN DEVELOPED BY A THIRD PARTY, THE THIRD PARTY
        !            45:  * DEVELOPER SHALL HAVE NO LIABILITY IN CONNECTION WITH THE USE,
        !            46:  * PERFORMANCE OR NON-PERFORMANCE OF THE SOFTWARE.
        !            47:  *
        !            48:  * $OpenXM: OpenXM_contrib2/asir2000/engine/mat.c,v 1.1.1.1 1999/12/03 07:39:08 noro Exp $
        !            49: */
1.1       noro       50: #include "ca.h"
                     51:
                     52: void addmat(vl,a,b,c)
                     53: VL vl;
                     54: MAT a,b,*c;
                     55: {
                     56:        int row,col,i,j;
                     57:        MAT t;
                     58:        pointer *ab,*bb,*tb;
                     59:
                     60:        if ( !a )
                     61:                *c = b;
                     62:        else if ( !b )
                     63:                *c = a;
                     64:        else if ( (a->row != b->row) || (a->col != b->col) ) {
                     65:                *c = 0; error("addmat : size mismatch");
                     66:        } else {
                     67:                row = a->row; col = a->col;
                     68:                MKMAT(t,row,col);
                     69:                for ( i = 0; i < row; i++ )
                     70:                        for ( j = 0, ab = BDY(a)[i], bb = BDY(b)[i], tb = BDY(t)[i];
                     71:                                j < col; j++ )
                     72:                                addr(vl,(Obj)ab[j],(Obj)bb[j],(Obj *)&tb[j]);
                     73:                *c = t;
                     74:        }
                     75: }
                     76:
                     77: void submat(vl,a,b,c)
                     78: VL vl;
                     79: MAT a,b,*c;
                     80: {
                     81:        int row,col,i,j;
                     82:        MAT t;
                     83:        pointer *ab,*bb,*tb;
                     84:
                     85:        if ( !a )
                     86:                chsgnmat(b,c);
                     87:        else if ( !b )
                     88:                *c = a;
                     89:        else if ( (a->row != b->row) || (a->col != b->col) ) {
                     90:                *c = 0; error("submat : size mismatch");
                     91:        } else {
                     92:                row = a->row; col = a->col;
                     93:                MKMAT(t,row,col);
                     94:                for ( i = 0; i < row; i++ )
                     95:                        for ( j = 0, ab = BDY(a)[i], bb = BDY(b)[i], tb = BDY(t)[i];
                     96:                                j < col; j++ )
                     97:                                subr(vl,(Obj)ab[j],(Obj)bb[j],(Obj *)&tb[j]);
                     98:                *c = t;
                     99:        }
                    100: }
                    101:
                    102: void mulmat(vl,a,b,c)
                    103: VL vl;
                    104: Obj a,b,*c;
                    105: {
                    106:        if ( !a || !b )
                    107:                *c = 0;
                    108:        else if ( OID(a) <= O_R )
                    109:                mulrmat(vl,(Obj)a,(MAT)b,(MAT *)c);
                    110:        else if ( OID(b) <= O_R )
                    111:                mulrmat(vl,(Obj)b,(MAT)a,(MAT *)c);
                    112:        else
                    113:                switch ( OID(a) ) {
                    114:                        case O_VECT:
                    115:                                switch ( OID(b) ) {
                    116:                                        case O_MAT:
                    117:                                                mulvectmat(vl,(VECT)a,(MAT)b,(VECT *)c); break;
                    118:                                        case O_VECT: default:
                    119:                                                notdef(vl,a,b,c); break;
                    120:                                }
                    121:                                break;
                    122:                        case O_MAT:
                    123:                                switch ( OID(b) ) {
                    124:                                        case O_VECT:
                    125:                                                mulmatvect(vl,(MAT)a,(VECT)b,(VECT *)c); break;
                    126:                                        case O_MAT:
                    127:                                                mulmatmat(vl,(MAT)a,(MAT)b,(MAT *)c); break;
                    128:                                        default:
                    129:                                                notdef(vl,a,b,c); break;
                    130:                                }
                    131:                                break;
                    132:                        default:
                    133:                                notdef(vl,a,b,c); break;
                    134:                }
                    135: }
                    136:
                    137: void divmat(vl,a,b,c)
                    138: VL vl;
                    139: Obj a,b,*c;
                    140: {
                    141:        Obj t;
                    142:
                    143:        if ( !b )
                    144:                error("divmat : division by 0");
                    145:        else if ( !a )
                    146:                *c = 0;
                    147:        else if ( OID(b) > O_R )
                    148:                notdef(vl,a,b,c);
                    149:        else {
                    150:                divr(vl,(Obj)ONE,b,&t); mulrmat(vl,t,(MAT)a,(MAT *)c);
                    151:        }
                    152: }
                    153:
                    154: void chsgnmat(a,b)
                    155: MAT a,*b;
                    156: {
                    157:        MAT t;
                    158:        int row,col,i,j;
                    159:        pointer *ab,*tb;
                    160:
                    161:        if ( !a )
                    162:                *b = 0;
                    163:        else {
                    164:                row = a->row; col = a->col;
                    165:                MKMAT(t,row,col);
                    166:                for ( i = 0; i < row; i++ )
                    167:                        for ( j = 0, ab = BDY(a)[i], tb = BDY(t)[i];
                    168:                                j < col; j++ )
                    169:                                chsgnr((Obj)ab[j],(Obj *)&tb[j]);
                    170:                *b = t;
                    171:        }
                    172: }
                    173:
                    174: void pwrmat(vl,a,r,c)
                    175: VL vl;
                    176: MAT a;
                    177: Obj r;
                    178: MAT *c;
                    179: {
                    180:        if ( !a )
                    181:                *c = 0;
                    182:        else if ( !r || !NUM(r) || !RATN(r) ||
                    183:                !INT(r) || (SGN((Q)r)<0) || (PL(NM((Q)r))>1) ) {
                    184:                *c = 0; error("pwrmat : invalid exponent");
                    185:        } else if ( a->row != a->col ) {
                    186:                *c = 0; error("pwrmat : non square matrix");
                    187:        }  else
                    188:                pwrmatmain(vl,a,QTOS((Q)r),c);
                    189: }
                    190:
                    191: void pwrmatmain(vl,a,e,c)
                    192: VL vl;
                    193: MAT a;
                    194: int e;
                    195: MAT *c;
                    196: {
                    197:        MAT t,s;
                    198:
                    199:        if ( e == 1 ) {
                    200:                *c = a;
                    201:                return;
                    202:        }
                    203:
                    204:        pwrmatmain(vl,a,e/2,&t);
                    205:        mulmat(vl,(Obj)t,(Obj)t,(Obj *)&s);
                    206:        if ( e % 2 )
                    207:                mulmat(vl,(Obj)s,(Obj)a,(Obj *)c);
                    208:        else
                    209:                *c = s;
                    210: }
                    211:
                    212: void mulrmat(vl,a,b,c)
                    213: VL vl;
                    214: Obj a;
                    215: MAT b,*c;
                    216: {
                    217:        int row,col,i,j;
                    218:        MAT t;
                    219:        pointer *bb,*tb;
                    220:
                    221:        if ( !a || !b )
                    222:                *c = 0;
                    223:        else {
                    224:                row = b->row; col = b->col;
                    225:                MKMAT(t,row,col);
                    226:                for ( i = 0; i < row; i++ )
                    227:                        for ( j = 0, bb = BDY(b)[i], tb = BDY(t)[i];
                    228:                                j < col; j++ )
                    229:                                mulr(vl,(Obj)a,(Obj)bb[j],(Obj *)&tb[j]);
                    230:                *c = t;
                    231:        }
                    232: }
                    233:
                    234: void mulmatmat(vl,a,b,c)
                    235: VL vl;
                    236: MAT a,b,*c;
                    237: {
                    238:        int arow,bcol,i,j,k,m;
                    239:        MAT t;
                    240:        pointer s,u,v;
                    241:        pointer *ab,*tb;
                    242:
                    243:        if ( a->col != b->row ) {
                    244:                *c = 0; error("mulmat : size mismatch");
                    245:        } else {
                    246:                arow = a->row; m = a->col; bcol = b->col;
                    247:                MKMAT(t,arow,bcol);
                    248:                for ( i = 0; i < arow; i++ )
                    249:                        for ( j = 0, ab = BDY(a)[i], tb = BDY(t)[i]; j < bcol; j++ ) {
                    250:                                for ( k = 0, s = 0; k < m; k++ ) {
                    251:                                        mulr(vl,(Obj)ab[k],(Obj)BDY(b)[k][j],(Obj *)&u); addr(vl,(Obj)s,(Obj)u,(Obj *)&v); s = v;
                    252:                                }
                    253:                                tb[j] = s;
                    254:                        }
                    255:                *c = t;
                    256:        }
                    257: }
                    258:
                    259: void mulmatvect(vl,a,b,c)
                    260: VL vl;
                    261: MAT a;
                    262: VECT b;
                    263: VECT *c;
                    264: {
                    265:        int arow,i,j,m;
                    266:        VECT t;
                    267:        pointer s,u,v;
                    268:        pointer *ab;
                    269:
                    270:        if ( !a || !b )
                    271:                *c = 0;
                    272:        else if ( a->col != b->len ) {
                    273:                *c = 0; error("mulmatvect : size mismatch");
                    274:        } else {
                    275:                arow = a->row; m = a->col;
                    276:                MKVECT(t,arow);
                    277:                for ( i = 0; i < arow; i++ ) {
                    278:                        for ( j = 0, s = 0, ab = BDY(a)[i]; j < m; j++ ) {
                    279:                                mulr(vl,(Obj)ab[j],(Obj)BDY(b)[j],(Obj *)&u); addr(vl,(Obj)s,(Obj)u,(Obj *)&v); s = v;
                    280:                        }
                    281:                        BDY(t)[i] = s;
                    282:                }
                    283:                *c = t;
                    284:        }
                    285: }
                    286:
                    287: void mulvectmat(vl,a,b,c)
                    288: VL vl;
                    289: VECT a;
                    290: MAT b;
                    291: VECT *c;
                    292: {
                    293:        int bcol,i,j,m;
                    294:        VECT t;
                    295:        pointer s,u,v;
                    296:
                    297:        if ( !a || !b )
                    298:                *c = 0;
                    299:        else if ( a->len != b->row ) {
                    300:                *c = 0; error("mulvectmat : size mismatch");
                    301:        } else {
                    302:                bcol = b->col; m = a->len;
                    303:                MKVECT(t,bcol);
                    304:                for ( j = 0; j < bcol; j++ ) {
                    305:                        for ( i = 0, s = 0; i < m; i++ ) {
                    306:                                mulr(vl,(Obj)BDY(a)[i],(Obj)BDY(b)[i][j],(Obj *)&u); addr(vl,(Obj)s,(Obj)u,(Obj *)&v); s = v;
                    307:                        }
                    308:                        BDY(t)[j] = s;
                    309:                }
                    310:                *c = t;
                    311:        }
                    312: }
                    313:
                    314: int compmat(vl,a,b)
                    315: VL vl;
                    316: MAT a,b;
                    317: {
                    318:        int i,j,t,row,col;
                    319:
                    320:        if ( !a )
                    321:                return b?-1:0;
                    322:        else if ( !b )
                    323:                return 1;
                    324:        else if ( a->row != b->row )
                    325:                return a->row>b->row ? 1 : -1;
                    326:        else if (a->col != b->col )
                    327:                return a->col > b->col ? 1 : -1;
                    328:        else {
                    329:                row = a->row; col = a->col;
                    330:                for ( i = 0; i < row; i++ )
                    331:                        for ( j = 0; j < col; j++ )
                    332:                                if ( t = compr(vl,(Obj)BDY(a)[i][j],(Obj)BDY(b)[i][j]) )
                    333:                                        return t;
                    334:                return 0;
                    335:        }
                    336: }
                    337:
                    338: pointer **almat_pointer(n,m)
                    339: int n,m;
                    340: {
                    341:        pointer **mat;
                    342:        int i;
                    343:
                    344:        mat = (pointer **)MALLOC(n*sizeof(pointer *));
                    345:        for ( i = 0; i < n; i++ )
                    346:                mat[i] = (pointer *)CALLOC(m,sizeof(pointer));
                    347:        return mat;
                    348: }

FreeBSD-CVSweb <freebsd-cvsweb@FreeBSD.org>