summaryrefslogtreecommitdiff
path: root/Trivac/src/NSSF.f
blob: 1470e220a367b848ccbde8fe372de1d8406eca92 (plain)
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
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
*DECK NSSF
      SUBROUTINE NSSF(NENTRY,HENTRY,IENTRY,JENTRY,KENTRY)
*
*-----------------------------------------------------------------------
*
*Purpose:
* Flux solution for a nodal method (NEM or ANM).
*
*Copyright:
* Copyright (C) 2021 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): A. Hebert
*
*Parameters: input/output
* NENTRY  number of LCM objects or files used by the operator.
* HENTRY  name of each LCM object or file:
*         HENTRY(1): create type(L_FLUX) nodal flux;
*         HENTRY(2): read-only type(L_TRACK) nodal tracking;
*         HENTRY(3): read-only type(L_MACROLIB) nodal macrolib.
* IENTRY  type of each LCM object or file:
*         =1 LCM memory object; =2 XSM file; =3 sequential binary file;
*         =4 sequential ascii file.
* JENTRY  access of each LCM object or file:
*         =0 the LCM object or file is created;
*         =1 the LCM object or file is open for modifications;
*         =2 the LCM object or file is open in read-only mode.
* KENTRY  LCM object address or file unit number.
*
*-----------------------------------------------------------------------
*
      USE GANLIB
*----
*  SUBROUTINE ARGUMENTS
*----
      INTEGER      NENTRY,IENTRY(NENTRY),JENTRY(NENTRY)
      TYPE(C_PTR)  KENTRY(NENTRY)
      CHARACTER    HENTRY(NENTRY)*12
*----
*  LOCAL VARIABLES
*----
      PARAMETER (NSTATE=40)
      TYPE(C_PTR) IPFLX,IPTRK,IPMAC
      CHARACTER HSIGN*12,TEXT4*4,TEXT12*12,HSMG*131,BNDTL*12
      LOGICAL LNODF
      INTEGER ISTATE(NSTATE)
      REAL REALIR
      DOUBLE PRECISION DBLLIR
      INTEGER, ALLOCATABLE, DIMENSION(:,:) :: ITRIAL
*----
*  PARAMETER VALIDATION.
*----
      IF(NENTRY.NE.3) CALL XABORT('NSSF: 3 PARAMETERS EXPECTED.')
      IF((IENTRY(1).NE.1).AND.(IENTRY(1).NE.2)) CALL XABORT('NSSF: LCM'
     1 //' OBJECT EXPECTED AT LHS.')
      DO IEN=2,3
        IF((IENTRY(IEN).NE.1).AND.(IENTRY(IEN).NE.2)) CALL XABORT('NSS'
     1  //'F: LCM OBJECT EXPECTED AT LHS.')
        IF(JENTRY(IEN).NE.2) CALL XABORT('NSSF: ENTRY IN READ-ONLY MOD'
     1  //'E EXPECTED.')
        CALL LCMGTC(KENTRY(IEN),'SIGNATURE',12,HSIGN)
        TEXT12=HENTRY(IEN)
        IF(IEN.EQ.2) THEN
          IF(HSIGN.NE.'L_TRACK') THEN
            CALL XABORT('NSSF: SIGNATURE OF '//TEXT12//' IS '//HSIGN//
     1      '. L_TRACK EXPECTED.')
          ENDIF
          IPTRK=KENTRY(2)
          CALL LCMPTC(KENTRY(1),'LINK.TRACK',12,HENTRY(2))
        ELSE IF(IEN.EQ.3) THEN
          IF(HSIGN.NE.'L_MACROLIB') THEN
            CALL XABORT('NSSF: SIGNATURE OF '//TEXT12//' IS '//HSIGN//
     1      '. L_MACROLIB EXPECTED.')
          ENDIF
          IPMAC=KENTRY(3)
          CALL LCMPTC(KENTRY(1),'LINK.MACRO',12,HENTRY(3))
        ENDIF
      ENDDO
      CALL LCMGTC(IPTRK,'TRACK-TYPE',12,TEXT12)
      IF(TEXT12.NE.'TRIVAC') CALL XABORT('NSSF: TRIVAC TRACKING EXPECT'
     1 //'ED.')
*----
*  PROCESS MACROLIB AND TRACKING.
*----
      CALL LCMGET(IPTRK,'STATE-VECTOR',ISTATE)
      NEL=ISTATE(1)
      NUN=ISTATE(2)
      NMIX=ISTATE(4)
      NADI=ISTATE(33)
      IGMAX=ISTATE(39)
      ICHX=ISTATE(12)
      IF((ICHX.LT.4).OR.(ICHX.GT.6)) THEN
        CALL XABORT('NSSF: CMFD, NEM OR ANM DISCRETIZATION EXPECTED.')
      ENDIF
      IF(ISTATE(6).EQ.2) THEN
        IDIM=1
      ELSE IF((ISTATE(6).EQ.5).AND.(ICHX.EQ.6)) THEN
        IDIM=2
      ELSE IF((ISTATE(6).EQ.7).AND.(ICHX.EQ.6)) THEN
        IDIM=3
      ELSE
        CALL XABORT('NSSF: 1D SLAB/2D-3D CARTESIAN GEOMETRY EXPECTED.')
      ENDIF
      IF(ISTATE(38).NE.0) CALL XABORT('NSSF: LUMP OPTION FORBIDDEN.')
      CALL LCMGET(IPMAC,'STATE-VECTOR',ISTATE)
      NG=ISTATE(1)
      IF(ISTATE(2).NE.NMIX) THEN
        WRITE(HSMG,'(39HNSSF: INVALID NUMBER OF MIXTURES (GEOM=,I5,
     1  10H MACROLIB=,I5,2H).)') NMIX,ISTATE(2)
        CALL XABORT(HSMG)
      ENDIF
*----
*  CREATE OR RECOVER THE FLUX.
*----
      IPFLX=KENTRY(1)
      IF(JENTRY(1).EQ.0) THEN
        HSIGN='L_FLUX'
        CALL LCMPTC(IPFLX,'SIGNATURE',12,HSIGN)
      ELSE IF(JENTRY(1).EQ.1) THEN
        CALL LCMGTC(IPFLX,'SIGNATURE',12,HSIGN)
        TEXT12=HENTRY(IEN)
        IF(HSIGN.NE.'L_FLUX') THEN
          CALL XABORT('NSSF: SIGNATURE OF '//TEXT12//' IS '//HSIGN//
     1    '. L_FLUX EXPECTED.')
        ENDIF
        CALL LCMGET(IPFLX,'STATE-VECTOR',ISTATE)
        IF(ISTATE(1).NE.NG) THEN
          WRITE(HSMG,'(41HNSSF: INVALID NUMBER OF GROUPS (MACROLIB=,I5,
     1    6H FLUX=,I5,2H).)') NG,ISTATE(1)
          CALL XABORT(HSMG)
        ELSE IF(ISTATE(2).NE.NUN) THEN
          WRITE(HSMG,'(43HNSSF: INVALID NUMBER OF UNKNOWNS (TRACKING=,
     1    I10,6H FLUX=,I10,2H).)') NUN,ISTATE(2)
          CALL XABORT(HSMG)
        ENDIF
      ELSE
        CALL XABORT('NSSF: FLUX IN CREATE OR MODIFICATION MODE EXPECTE'
     1  //'D.')
      ENDIF
*---
*  READ DATA
*---
      ALLOCATE(ITRIAL(NMIX,NG))
      IPRINT=1
      ICL1=3
      ICL2=3
      MAXNOD=300
      MAXTHR=0
      MAXOUT=100
      EPSNOD=1.0E-6
      EPSTHR=1.0E-6
      EPSOUT=1.0E-5
      LNODF=.FALSE.
      BB2=0.0
      BNDTL='quadratic'
      NPASS=3
      ITRIAL(:,:)=1
   10 CALL REDGET(ITYPLU,INTLIR,REALIR,TEXT4,DBLLIR)
      IF(ITYPLU.EQ.10) GO TO 100
   20 IF(ITYPLU.NE.3) CALL XABORT('NSSF: READ ERROR - CHARACTER VARIAB'
     > //'LE EXPECTED')
      IF(TEXT4.EQ.';') THEN
        GO TO 100
      ELSE IF(TEXT4.EQ.'EDIT') THEN
        CALL REDGET(ITYPLU,IPRINT,REALIR,TEXT4,DBLLIR)
        IF(ITYPLU.NE.1) CALL XABORT('NSSF: INTEGER DATA EXPECTED(1).')
      ELSE IF((TEXT4.EQ.'VAR1').OR.(TEXT4.EQ.'ACCE')) THEN
        CALL REDGET(ITYPLU,ICL1,REALIR,TEXT4,DBLLIR)
        IF(ITYPLU.NE.1) CALL XABORT('NSSF: INTEGER DATA EXPECTED(2).')
        CALL REDGET(ITYPLU,ICL2,REALIR,TEXT4,DBLLIR)
        IF(ITYPLU.NE.1) CALL XABORT('NSSF: INTEGER DATA EXPECTED(3).')
      ELSE IF(TEXT4=='ADI') THEN
        CALL REDGET(ITYPLU,NADI,REALIR,TEXT4,DBLLIR)
        IF(ITYPLU.NE.1) CALL XABORT('NSSF: INTEGER DATA EXPECTED(5).')
      ELSE IF(TEXT4=='NUPD') THEN
        ! maximum number of nodal correction iterations
   30   CALL REDGET(ITYPLU,INTLIR,REALIR,TEXT4,DBLLIR)
        IF(ITYPLU.EQ.1) THEN
          MAXNOD=INTLIR
          CALL REDGET(ITYPLU,INTLIR,REALIR,TEXT4,DBLLIR)
          IF(ITYPLU.EQ.1) THEN
            NPASS=INTLIR
            GO TO 30
          ENDIF
        ELSE IF(ITYPLU.EQ.2) THEN
          EPSNOD=REALIR
        ELSE
          GO TO 20
        ENDIF
        GO TO 30
      ELSE IF(TEXT4=='EXTE') THEN
        ! maximum number and convergence criterion of Keff iterations
   40   CALL REDGET(ITYPLU,INTLIR,REALIR,TEXT4,DBLLIR)
        IF(ITYPLU.EQ.1) THEN
          MAXOUT=INTLIR
        ELSE IF(ITYPLU.EQ.2) THEN
          EPSOUT=REALIR
        ELSE
          GO TO 20
        ENDIF
        GO TO 40
      ELSE IF(TEXT4=='THER') THEN
        ! maximum number and convergence criterion of thermal iterations
   50   CALL REDGET(ITYPLU,INTLIR,REALIR,TEXT4,DBLLIR)
        IF(ITYPLU.EQ.1) THEN
          MAXTHR=INTLIR
        ELSE IF(ITYPLU.EQ.2) THEN
          EPSTHR=REALIR
        ELSE
          GO TO 20
        ENDIF
        GO TO 50
      ELSE IF(TEXT4.EQ.'NODF') THEN
        LNODF=.TRUE.
      ELSE IF(TEXT4=='LEAK') THEN
        CALL REDGET(ITYPLU,INTLIR,REALIR,BNDTL,DBLLIR)
        IF(ITYPLU/=3) CALL XABORT('NSSF: READ ERROR - CHARACTER VARIAB'
     >  //'LE EXPECTED')
        IF((BNDTL.NE.'flat').AND.(BNDTL.NE.'quadratic')) THEN
         CALL XABORT('NSSF: flat OR quadratic KEYWORD EXPECTED')
        ENDIF
      ELSE IF(TEXT4.EQ.'BUCK') THEN
        CALL REDGET(ITYPLU,INTLIR,BB2,TEXT4,DBLLIR)
        IF(ITYPLU.NE.2) CALL XABORT('NSSF: READ ERROR - REAL VARIABLE '
     >  //'EXPECTED')
      ELSE
        CALL XABORT('NSSF: ILLEGAL KEYWORD '//TEXT4)
      ENDIF
      GO TO 10
  100 IF(IGMAX.GT.NG) CALL XABORT('NSSF: IGMAX>NG.')
      IF(IPRINT.GT.0) THEN
        WRITE(6,'(/47H NSSF: number of transverse current iterations=,
     >  I3)') NPASS
      ENDIF
      IF(IGMAX.GT.0) ITRIAL(:NMIX,IGMAX:NG)=2
      CALL NSSDRV(IPTRK,IPMAC,IPFLX,ICHX,IDIM,NUN,NG,NEL,NMIX,ITRIAL,
     > ICL1,ICL2,NADI,EPSNOD,MAXNOD,EPSTHR,MAXTHR,EPSOUT,MAXOUT,LNODF,
     > BNDTL,NPASS,BB2,IPRINT)
      DEALLOCATE(ITRIAL)
      RETURN
      END