1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
|
*DECK MTLDLF
SUBROUTINE MTLDLF(NAMP,IPTRK,IPSYS,ITY,IMPX)
*
*-----------------------------------------------------------------------
*
*Purpose:
* LCM driver for the L-D-L(t) factorization of a symmetric matrix.
* The factorized matrix is stored on LCM under name 'I'//NAMP.
*
*Copyright:
* Copyright (C) 2002 Ecole Polytechnique de Montreal
* This library is free software; you can redistribute it and/or
* modify it under the terms of the GNU Lesser General Public
* License as published by the Free Software Foundation; either
* version 2.1 of the License, or (at your option) any later version
*
*Author(s): Alain Hebert
*
*Parameters: input
* NAMP name of the coefficient matrix.
* IPTRK L_TRACK pointer to the tracking information.
* IPSYS L_SYSTEM pointer to system matrices.
* ITY type of coefficient matrix (1: Bivac; 2: classical Trivac;
* 3: Thomas-Raviart).
* IMPX print flag (equal to zero for no print).
*
*-----------------------------------------------------------------------
*
USE GANLIB
*----
* SUBROUTINE ARGUMENTS
*----
TYPE(C_PTR) IPTRK,IPSYS
CHARACTER NAMP*12
INTEGER ITY,IMPX
*----
* LOCAL VARIABLES
*----
PARAMETER (NSTATE=40,NPREF=5)
CHARACTER HIN*12,HOUT*12,PREFIX(5)*2,NAMLCM*12,NAMMY*12
LOGICAL EMPTY,LCM
INTEGER ITP(NSTATE)
INTEGER, DIMENSION(:), ALLOCATABLE :: MU,NBL,LBL
REAL, DIMENSION(:), ALLOCATABLE :: T
REAL, DIMENSION(:), POINTER :: ASM
TYPE(C_PTR) ASM_PTR
DATA (PREFIX(I),I=1,NPREF)/' ','W_','X_','Y_','Z_'/
*
IF(ITY.EQ.1) THEN
* BIVAC TRACKING.
CALL LCMGET(IPTRK,'STATE-VECTOR',ITP)
ISEG=0
NLF=ITP(14)
ELSE
* CLASSICAL TRIVAC TRACKING.
CALL LCMGET(IPTRK,'STATE-VECTOR',ITP)
ISEG=ITP(17)
NLF=ITP(30)
ENDIF
*
DO 30 IS=1,NPREF
IF(PREFIX(IS).EQ.' ') THEN
HIN=NAMP
HOUT='I'//NAMP(:11)
ELSE
HIN=PREFIX(IS)//NAMP(:11)
HOUT=PREFIX(IS)(:1)//'I'//NAMP(:10)
ENDIF
*----
* PERFORM FACTORIZATION OF MATRICES
*----
CALL LCMLEN(IPSYS,HIN,ILENG,ITYLCM)
IF(ILENG.GT.0) THEN
IF(ISEG.EQ.0) THEN
CALL LCMLEN(IPTRK,'MU'//PREFIX(IS)(:1),LMU,ITYLCM)
ALLOCATE(MU(LMU))
CALL LCMGET(IPTRK,'MU'//PREFIX(IS)(:1),MU)
ELSE
CALL LCMLEN(IPTRK,'MUV'//PREFIX(IS)(:1),LMU,ITYLCM)
ALLOCATE(MU(LMU))
CALL LCMGET(IPTRK,'MUV'//PREFIX(IS)(:1),MU)
ENDIF
ILEN=MU(LMU)
IF(NLF.GT.0) ILEN=ILEN*NLF/2
IF(IMPX.GT.0) THEN
CALL LCMINF(IPSYS,NAMLCM,NAMMY,EMPTY,ILONG,LCM)
WRITE(6,'(/30H MTLDLF: FACTORIZATION OF LCM ,
1 8HMATRIX '',A12,23H''. CREATION OF MATRIX '',A12,
2 14H'' LOCATED IN '',A12,2H''.)') HIN,HOUT,NAMLCM
ENDIF
ASM_PTR=LCMARA(ILENG)
CALL C_F_POINTER(ASM_PTR,ASM,(/ ILENG /))
CALL LCMGET(IPSYS,HIN,ASM)
IF(ISEG.EQ.0) THEN
IF(ILEN.NE.ILENG) CALL XABORT('MTLDLF: INCONSISTENT INF'
1 //'ORMATION ON LCM (1).')
IF(NLF.EQ.0) THEN
CALL ALLDLF(LMU,ASM(1),MU)
ELSE
IOF=1
DO 10 IL=0,NLF-2,2
CALL ALLDLF(LMU,ASM(IOF),MU)
IOF=IOF+MU(LMU)
10 CONTINUE
ENDIF
ELSE
IF(ISEG*ILEN.NE.ILENG) CALL XABORT('MTLDLF: INCONSISTEN'
1 //'T INFORMATION ON LCM (2).')
CALL LCMLEN(IPTRK,'NBL'//PREFIX(IS)(:1),LON,ITYLCM)
ALLOCATE(NBL(LON),LBL(LON))
CALL LCMGET(IPTRK,'NBL'//PREFIX(IS)(:1),NBL)
CALL LCMGET(IPTRK,'LBL'//PREFIX(IS)(:1),LBL)
ALLOCATE(T(ISEG))
IF(NLF.EQ.0) THEN
CALL ALVDLF(ASM(1),MU,ISEG,LON,NBL,LBL,T)
ELSE
IOF=1
DO 20 IL=0,NLF-2,2
CALL ALVDLF(ASM(IOF),MU,ISEG,LON,NBL,LBL,T)
IOF=IOF+MU(LMU)
20 CONTINUE
ENDIF
DEALLOCATE(T,LBL,NBL)
ENDIF
DEALLOCATE(MU)
CALL LCMPPD(IPSYS,HOUT,ILENG,2,ASM_PTR)
ENDIF
30 CONTINUE
RETURN
END
|