summaryrefslogtreecommitdiff
path: root/Donjon/src/RODMOD.f
blob: 25cd8e92d209db00a2cfc1b555aed13d251c36d8 (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
*DECK RODMOD
      SUBROUTINE RODMOD(IPMAP,NGRP,MAXMIX,NCH,NB,RNAME,INS,INSS,
     > HMIX,RMIX,NREB,RODSIZE,INTROD,INI,NUMMIX,NCALL)
*
*-----------------------------------------------------------------------
*
*Purpose:
* Modification of the rod data stored in the PARAM folder of a fuel map
*
*Copyright:
* Copyright (C) 2017 Ecole Polytechnique de Montreal.
*
*Author(s): 
* G. Tixier
*
*Parameters: input
* IPMAP   pointer to the fuel map
* NGRP    number of rod groups
* MAXMIX  maximum number of rod mix
* NCH     number of fuel channels
* NB      number of fuel bundles per channel
* RNAME   name of rod group
* INS     rod insertion for each rod group
* INSS    rod insertion step (in cm)
* HMIX    height of rod mix (if more than 2 rod mix are defined)
* RMIX    number associated to rod mix
* NREB    number of bottom-reflective meshes
* RODSIZE total number of meshes for the fuel zone (=NCH*NB)
* INTROD  final value of data stored in the PARAM folder after 
*         calculation
* INI     initial value of data stored in the PARAM folder (no rod
*         inserted)
* NUMMIX  number of rod mix for each rod group
* NCALL   number to distinguish first or other calls to the ROD: 
*         module
*
*-----------------------------------------------------------------------
*
      USE GANLIB
*----
*  SUBROUTINE ARGUMENTS
*----
      INTEGER NGRP,MAXMIX,NREB,NCH,NB,RODSIZE,NUMMIX(NGRP)
      INTEGER INSM,NCALL
      REAL HMIX(NGRP*MAXMIX),RMIX(NGRP*MAXMIX),INTROD(RODSIZE),INS(NGRP)
      REAL INI,INSS
      CHARACTER(LEN=3) RNAME(NGRP),RASS(NCH)
      TYPE(C_PTR) IPMAP
*----
*  LOCAL VARIABLES
*----
      CHARACTER TEXT*3
      INTEGER I,J,K,L,M
      REAL FLOT,F1,F2,ICM,HMB,HMT,ZMI1,ZMI2
      DOUBLE PRECISION DFLOT
      TYPE(C_PTR) JPMAP,MPMAP
*----
*  ALLOCATABLE ARRAYS
*----
      REAL, ALLOCATABLE, DIMENSION(:) :: GEOZZ

* Recover axial meshing
      JPMAP=LCMGID(IPMAP,'GEOMAP')
      CALL LCMLEN(JPMAP,'MESHZ',ILONG,ITYLCM)
      ALLOCATE(GEOZZ(ILONG))
      CALL LCMGET(JPMAP,'MESHZ',GEOZZ)

* Recover rod parameters
      MPMAP=LCMGID(IPMAP,'ROD-INFO')
      CALL LCMGET(MPMAP,'ROD-INIT',INI)
      CALL LCMGET(MPMAP,'INS-MAX',INSM)
      CALL LCMGET(MPMAP,'STEP-CM',INSS)
      CALL LCMGET(MPMAP,'REFL-BOTTOM',NREB)      
      CALL LCMGTC(MPMAP,'ROD-NAME',3,NGRP,RNAME)
      CALL LCMGET(MPMAP,'ROD-INSERT',INS)
      CALL LCMGET(MPMAP,'ROD-RIN',RMIX)
      CALL LCMGET(MPMAP,'ROD-NBZONE',NUMMIX)
      CALL LCMGET(MPMAP,'ROD-HEIGHT',HMIX)
      IF(NB+NREB+1.GT.ILONG) CALL XABORT('RODMOD: GEOZZ OVERFLOW.')

      IF(NCALL.EQ.1) THEN
* If it is the first call to the ROD: module, the rod map is stored in
* the fuel map
      N=1
      CALL REDGET(ITYP,NITMA,FLOT,TEXT,DFLOT)
      IF(ITYP.NE.3)CALL XABORT('@RODMOD: CHARACTER DATA EXPECTED.')
      DO WHILE(N.LE.NCH)
        RASS(N)=TEXT
        CALL REDGET(ITYP,NITMA,FLOT,TEXT,DFLOT)
        N=N+1
      END DO
      CALL LCMPTC(MPMAP,'ROD-MAP',3,NCH,RASS)
      ENDIF
* Recover rod map
      CALL LCMGTC(MPMAP,'ROD-MAP',3,NCH,RASS)
      I=1
      DO WHILE(I.LE.NCH)
        K=NB
        DO WHILE(K.GE.1)
        IF((RASS(I).EQ.'|').OR.(RASS(I).EQ.'-').OR.(RASS(I).EQ.'-|-'
     >  )) THEN
* If no control rod is defined
          INTROD((K-1)*NCH+I)=INI
          GO TO 10
        ELSE
* Recover control rod information
          J=1
          DO WHILE(J.LE.NGRP)
            IF(RASS(I).EQ.RNAME(J)) THEN
              EXIT
            ELSE
              J=J+1
            ENDIF 
          END DO
          M=NUMMIX(J)
          DO WHILE(M.GE.1)
            IF(M.EQ.1) THEN
              ICM=INS(J)*INSS
            ELSE
              ICM=INS(J)*INSS-HMIX(J+(M-2)*NGRP)
            ENDIF
* Mesh size calculations
            HMB=GEOZZ(NB+NREB+1)-GEOZZ(K+NREB)
            HMT=GEOZZ(NB+NREB+1)-GEOZZ(K+NREB+1)
            IF(ICM.LT.0.0) THEN
              M=M-1
            ELSE
              ZMI1=HMB-ICM
              IF(ZMI1.LE.0.0) THEN
                INTROD((K-1)*NCH+I)=RMIX(J+(M-1)*NGRP)
                EXIT
              ELSE
                ZMI2=HMT-ICM
                IF(ZMI2.GE.0.0) THEN
                  IF(M.EQ.1) THEN
                    INTROD((K-1)*NCH+I)=INI
                    GO TO 10
                  ELSE
                    M=M-1
                  ENDIF
                ELSE
* Calculation of the proportion of control rod for the mesh considered
                    F1=HMB-ICM
                    F2=GEOZZ(K+NREB+1)-GEOZZ(K+NREB)-F1
                  IF(M.EQ.1) THEN
                    INTROD((K-1)*NCH+I)=(F2*RMIX(J+(M-1)*NGRP)
     >              +F1*INI)/(F1+F2)
                    GO TO 10
                  ELSE
                    IF(M.GE.3) THEN
* It is not possible to have two interfaces in one mesh
                    IF(((INS(J)*INSS-HMIX(J+(M-3)*NGRP)).GE.HMT).AND.
     >              ((INS(J)*INSS-HMIX(J+(M-3)*NGRP)).LE.HMB))
     >              CALL XABORT('@RODMOD: IT IS NOT POSSIBLE TO HAVE'
     1   //' A MESH WITH MORE THAN TWO INTERFACES. HMIX HAS TO BE'
     1   //' BIGGER')
                    ENDIF
* If two mixtures fill one mesh, we have to compute a fraction of
* insertion for each rod mixture
                    INTROD((K-1)*NCH+I)=(F2*RMIX(J+(M-1)*NGRP)+
     >              F1*RMIX(J+(M-2)*NGRP))/(F1+F2)
                    EXIT
                  ENDIF
                ENDIF
              ENDIF
            ENDIF
          END DO
        ENDIF
        K=K-1
        END DO
   10   L=K-1
* If a control rod does not fill one mesh entirely, the meshes located
* at the bottom of this mesh do not contain control rod
        DO WHILE(L.GE.1)
          INTROD((L-1)*NCH+I)=INI
          L=L-1
        END DO
      I=I+1
      END DO
      DEALLOCATE(GEOZZ)
      RETURN
      END