MINICERN/packlib/kernlib/kernnum/f010fort/cinv.F

   1 *
   2 * $Id$
   3 *
   4 * $Log$
   5 * Revision 1.1.1.1  1996/02/15 17:48:49  mclareni
   6 * Kernlib
   7 *
   8 *
   9 #include "kernnum/pilot.h"
  10       SUBROUTINE CINV(N,A,IDIM,R,IFAIL)
  11       REAL R(N),T1,T2,T3
  12       COMPLEX A(IDIM,N),ONE,DET,TEMP,S,
  13      $        C11,C12,C13,C21,C22,C23,C31,C32,C33
  14       CHARACTER*6 NAME
  15       DATA NAME/'CINV'/,KPRNT/0/
  16       DATA ONE/(1.0,0.0)/
  17 C
  18 C     ******************************************************************
  19 C
  20 C     REPLACES A BY ITS INVERSE.
  21 C
  22 C     (PARAMETERS AS FOR CEQINV.)
  23 C
  24 C     CALLS ... CFACT, CFINV, F010PR, ABEND.
  25 C
  26 C     ******************************************************************
  27 C
  28 C  TEST FOR PARAMETER ERRORS.
  29 C
  30       IF((N.LT.1).OR.(N.GT.IDIM)) GO TO 7
  31 C
  32 C  TEST FOR N.LE.3.
  33 C
  34       IF(N.GT.3) GO TO 6
  35       IFAIL=0
  36       IF(N.LT.3) GO TO 4
  37 C
  38 C  N=3 CASE.
  39 C
  40 C     COMPUTE COFACTORS.
  41       C11=A(2,2)*A(3,3)-A(2,3)*A(3,2)
  42       C12=A(2,3)*A(3,1)-A(2,1)*A(3,3)
  43       C13=A(2,1)*A(3,2)-A(2,2)*A(3,1)
  44       C21=A(3,2)*A(1,3)-A(3,3)*A(1,2)
  45       C22=A(3,3)*A(1,1)-A(3,1)*A(1,3)
  46       C23=A(3,1)*A(1,2)-A(3,2)*A(1,1)
  47       C31=A(1,2)*A(2,3)-A(1,3)*A(2,2)
  48       C32=A(1,3)*A(2,1)-A(1,1)*A(2,3)
  49       C33=A(1,1)*A(2,2)-A(1,2)*A(2,1)
  50       T1=ABS(REAL(A(1,1)))+ABS(AIMAG(A(1,1)))
  51       T2=ABS(REAL(A(2,1)))+ABS(AIMAG(A(2,1)))
  52       T3=ABS(REAL(A(3,1)))+ABS(AIMAG(A(3,1)))
  53 C
  54 C     (SET TEMP=PIVOT AND DET=PIVOT*DET.)
  55       IF(T1.GE.T2) GO TO 1
  56          IF(T3.GE.T2) GO TO 2
  57 C        (PIVOT IS A21)
  58             TEMP=A(2,1)
  59             DET=C13*C32-C12*C33
  60             GO TO 3
  61     1 IF(T3.GE.T1) GO TO 2
  62 C     (PIVOT IS A11)
  63          TEMP=A(1,1)
  64          DET=C22*C33-C23*C32
  65          GO TO 3
  66 C     (PIVOT IS A31)
  67     2    TEMP=A(3,1)
  68          DET=C23*C12-C22*C13
  69 C
  70 C     SET ELEMENTS OF INVERSE IN A.
  71     3 IF( REAL(DET).EQ.0. .AND. AIMAG(DET).EQ.0. ) GO TO 8
  72       S=TEMP/DET
  73       A(1,1)=S*C11
  74       A(1,2)=S*C21
  75       A(1,3)=S*C31
  76       A(2,1)=S*C12
  77       A(2,2)=S*C22
  78       A(2,3)=S*C32
  79       A(3,1)=S*C13
  80       A(3,2)=S*C23
  81       A(3,3)=S*C33
  82       RETURN
  83 C
  84     4 IF(N.LT.2) GO TO 5
  85 C
  86 C  N=2 CASE BY CRAMERS RULE.
  87 C
  88       DET=A(1,1)*A(2,2)-A(1,2)*A(2,1)
  89       IF( REAL(DET).EQ.0. .AND. AIMAG(DET).EQ.0. ) GO TO 8
  90       S=ONE/DET
  91       C11   =S*A(2,2)
  92       A(1,2)=-S*A(1,2)
  93       A(2,1)=-S*A(2,1)
  94       A(2,2)=S*A(1,1)
  95       A(1,1)=C11
  96       RETURN
  97 C
  98 C  N=1 CASE.
  99 C
 100     5 IF( REAL(A(1,1)).EQ.0. .AND. AIMAG(A(1,1)).EQ.0. ) GO TO 8
 101       A(1,1)=ONE/A(1,1)
 102       RETURN
 103 C
 104 C  N.GT.3 CASES.  FACTORIZE MATRIX AND INVERT.
 105 C
 106     6 CALL CFACT(N,A,IDIM,R,IFAIL,DET,JFAIL)
 107       IF(IFAIL.NE.0) RETURN
 108       CALL CFINV(N,A,IDIM,R)
 109       RETURN
 110 C
 111 C  ERROR EXITS.
 112 C
 113     7 IFAIL=+1
 114       CALL F010PR(NAME,N,IDIM,K,KPRNT)
 115       RETURN
 116 C
 117     8 IFAIL=-1
 118       RETURN
 119 C
 120       END