OpenRadioss 2025.1.11
OpenRadioss project
Loading...
Searching...
No Matches
freimpl.F File Reference
#include "implicit_f.inc"
#include "units_c.inc"
#include "impl1_c.inc"
#include "impl2_c.inc"
#include "parit_c.inc"
#include "com01_c.inc"
#include "com06_c.inc"
#include "buckcom.inc"
#include "scr06_c.inc"

Go to the source code of this file.

Functions/Subroutines

subroutine freimpl (ikad, key0, kimpl)
subroutine order_dtf (n, rc)

Function/Subroutine Documentation

◆ freimpl()

subroutine freimpl ( integer, dimension(0:*) ikad,
character, dimension(*) key0,
integer kimpl )

Definition at line 40 of file freimpl.F.

41C-----------------------------------------------
42C M o d u l e s
43C-----------------------------------------------
44 USE imp_dyna
45 USE imp_kbcs
46 USE imp_pcg_proj
47 USE imp_spbrm
48 USE message_mod
50C-----------------------------------------------
51C I m p l i c i t T y p e s
52C-----------------------------------------------
53#include "implicit_f.inc"
54C-----------------------------------------------
55C D u m m y A r g u m e n t s
56C-----------------------------------------------
57 INTEGER IKAD(0:*),KIMPL
58 CHARACTER KEY0(*)*5
59C-----------------------------------------------
60C C o m m o n B l o c k s
61C-----------------------------------------------
62#include "units_c.inc"
63#include "impl1_c.inc"
64#include "impl2_c.inc"
65#include "parit_c.inc"
66#include "com01_c.inc"
67#include "com06_c.inc"
68#include "buckcom.inc"
69#include "scr06_c.inc"
70C-----------------------------------------------
71C E x t e r n a l F u n c t i o n s
72C-----------------------------------------------
73 INTEGER NVAR
74C-----------------------------------------------
75C L o c a l V a r i a b l e s
76C-----------------------------------------------
77 INTEGER I, NBC, K, IKEY,IM,J,NJ,KK
78 CHARACTER TITLE*72, KEY2*5, KEY3*5, KEY4*5
79 CHARACTER(LEN=NCHARLINE100)::CARTE
80C----------------------------------------
81 ikey=kimpl
82 impl_s=0
83 idyna=0
84 iline=0
85 isprb=0
86 isolv=0
87 insolv=0
88 idtc=0
89 im=0
90 ikg=1
91 kz_tol=zero
92 sk_int=zero
93 d_tol=zero
94 lprint=0
95 nprint=0
96 impdeb=0
97 solvnfo=0
98 prstifmat=0
99 prstifmat_tol=zero
100 prstifmat_nc=1
101 prstifmat_it=0
102 impmv=1
103 isigini=0
104 ilintf=0
105 iprec = 0
106 l_lim = 0
107 itol = 0
108 l_tol =zero
109 dt_imp = zero
110 dt_min = zero
111 dt_max = zero
112 imp_rby=0
113 imp_int=0
114 isprn = 1
115C INTP_C = 0
116C -----after debugging on int24 spmd, change defaut to INTP_C=1 (INTP_C=0 suppressed good for maintenance)
117 intp_c = 1
118 l_bfgs = 0
119C IRREF = 0
120 irref = 1
121 iqstat = 0
122 ibuckl = 0
123 iscau = 0
124 imp_lr=0
125 ikproj=0
126 ismdisp = 0
127 IF(ikad(ikey)/=ikad(ikey+1))THEN
128 k=0
129 impl_s=1
130 ncinp=1
131 n_pat = 1
132 imp_chk = 0
133 imp_int7 = 0
134 ittoff = 0
135 scal_dtq = one
136 idy_damp=0
137 iautspc = 1
138 itrmax = 0
139 msg_lvl = 0
140 b_order =0
141 b_mcore =0
142 irefi = 0
143 iline_s = 0
144 nls_lim = 0
145 ls_tol = zero
146 ndiver = 0
147 ikt = 0
148 ndtfix = 0
149 ikpres = 1
150 n_tolu=zero
151 n_tolf=zero
152 n_tole=zero
153 ncy_max = 0
154 rf_min = zero
155 rf_max = zero
156 ipupd = 0
157 tol_div = zero
158 m_vs = 0
159 ipro_s0=0
160 iikgoff = 1
161 m_msg = 0
162 m_order =0
163 m_ocore =0
164 irig_m = 0
165 1160 READ(iusc1,rec=ikad(ikey)+k,fmt='(7X,A,1X,A,1X,A,25X,I10)',err=9990)key2,key3,key4,nbc
166 k=k+1
167C----------------------------
168C Dynamic implicit
169C----------------------------
170 IF(key2(1:4)=='DYNA')THEN
171 IF (idyna==0) idyna=1
172 IF(key3(1:4)=='DAMP')THEN
173 idy_damp=1
174 CALL wriusc2(ikad(ikey)+k,1,key0(ikey))
175 READ(iusc2,*) dampa_imp,dampb_imp
176 ELSE IF(key3(1:3)=='FSI')THEN
177 WRITE(6,*) "ERROR: /IMPL/DYNA/FSI IS A DEPRECATED FEATURE"
178 GOTO 9990
179 ELSE
180 READ(key3,'(I2)')im
181 idyna=max(idyna,im)
182 IF(idyna==1)THEN
183 CALL wriusc2(ikad(ikey)+k,1,key0(ikey))
184 READ(iusc2,*)hht_a
185 ELSEIF(idyna==2)THEN
186 CALL wriusc2(ikad(ikey)+k,1,key0(ikey))
187 READ(iusc2,*)newm_a,newm_b
188 ELSE
189 hht_a =-em20
190 ENDIF
191 ENDIF
192C----------------------------
193C Implicit linear
194C----------------------------
195 ELSEIF(key2(1:4)=='LINE')THEN
196 iline=1
197 IF(key3(1:5)=='INTER') THEN
198 READ(key4,'(I5)')ilintf
199 ilintf = max(2,ilintf)
200 ELSEIF(key3(1:5)=='SCAUC') THEN
201 iscau = 1
202 ENDIF
203 ELSEIF(key2(1:5)=='MONVO')THEN
204 IF(key3(1:3)=='OFF')impmv=0
205 ELSEIF(key2(1:5)=='SPRIN')THEN
206 IF(key3(1:4)=='NONL')THEN
207 isprn = 1
208 ELSEIF(key3(1:4)=='LINE')THEN
209 isprn = 0
210 ELSE
211 GOTO 9990
212 ENDIF
213 ELSEIF(key2(1:5)=='PREPA')THEN
214 READ(key3,'(I2)')n_pat
215 ELSEIF(key2(1:5)=='PROJV')THEN
216 READ(key3,'(I2)') m_vs
217 ELSEIF(key2(1:5)=='PROSI')THEN
218 READ(key3,'(I2)') ipro_s0
219C----------------------------
220C Implicit check
221C----------------------------
222 ELSEIF(key2(1:5)=='CHECK')THEN
223 imp_chk = 1
224C----------------------------
225C Implicit quasi-static
226C----------------------------
227 ELSEIF(key2(1:5)=='QSTAT')THEN
228 iqstat = 1
229 IF(key3(1:5)=='DTSCA')THEN
230 CALL wriusc2(ikad(ikey)+k,1,key0(ikey))
231 READ(iusc2,*)scal_dtq
232 ELSEIF(key3(1:5)=='MRIGM')THEN
233 irig_m = 1
234 CALL wriusc2(ikad(ikey)+k,1,key0(ikey))
235 READ(iusc2,*,err=520,END=520)E_REF(1),E_REF(2),E_REF(3)
236 IF (e_ref(1)>0.AND.e_ref(2)>0.AND.e_ref(3)>0) irig_m = 2
237 520 CONTINUE
238 ELSE
239 READ(key3,'(I2)')im
240 iqstat=max(iqstat,im)
241 ENDIF
242C----------------------------
243C spring-back
244C----------------------------
245 ELSEIF(key2(1:4)=='SPRB')THEN
246 isprb=1
247C----------------------------
248C print-out
249C----------------------------
250 ELSEIF(key2=='PRINT')THEN
251 IF(key3(1:4)=='LINE')THEN
252 READ(key4,'(I5)')lprint
253 ELSEIF(key3(1:4)=='NONL')THEN
254 READ(key4,'(I5)')nprint
255 ELSEIF(key3(1:4)=='STIF')THEN
256 prstifmat = 1
257 CALL wriusc2(ikad(ikey)+k,1,key0(ikey))
258 READ(iusc2,*)prstifmat_tol,prstifmat_nc,prstifmat_it
259 ELSE
260 GOTO 9990
261 ENDIF
262C----------------------------
263C Linear SOLVER
264C----------------------------
265 ELSEIF(key2(1:4)=='SOLV')THEN
266 READ(key3,'(I2)')isolv
267 CALL wriusc2(ikad(ikey)+k,1,key0(ikey))
268 READ(iusc2,*)iprec,l_lim,itol,l_tol
269 IF (isolv==3) imumpsd=l_lim
270C----------------------------
271C BCS parameters
272C----------------------------
273 ELSEIF(key2(1:4)=='SBCS')THEN
274 IF(key3(1:5)=='MSGLV')THEN
275 READ(key4,'(i2)')MSG_LVL
276 ELSEIF(KEY3(1:5)=='order')THEN
277 READ(KEY4,'(i2)')B_ORDER
278C-------0 default 1 :MMD 2 :metis
279 ELSEIF(KEY3(1:5)=='outco')THEN
280 B_MCORE=1
281 ELSE
282 GOTO 9990
283 ENDIF
284C----------------------------
285C MUMPS parameters
286C----------------------------
287 ELSEIF(KEY2(1:5)=='mumps')THEN
288 IF(KEY3(1:5)=='msglv')THEN
289 READ(KEY4,'(i2)')M_MSG
290 ELSEIF(KEY3(1:5)=='order')THEN
291 IF(KEY4(1:5)=='metis')THEN
292 M_ORDER = 5
293 ELSEIF(KEY4(1:4)=='pord')THEN
294 M_ORDER = 4
295 END IF
296C-------0 default 1 :MMD 2 :metis
297 ELSEIF(KEY3(1:5)=='outco')THEN
298 M_OCORE=1
299 ELSEIF(KEY3(1:5)=='autoc')THEN
300 M_OCORE=-1
301 ELSE
302 GOTO 9990
303 ENDIF
304C----------------------------
305C Nonlinear SOLVER
306C----------------------------
307 ELSEIF(KEY2(1:4)=='nonl')THEN
308 IF(KEY3(1:5)=='ktang')THEN
309 IKT = 1
310 ELSEIF(KEY3(1:5)=='ktful')THEN
311 IKT = 2
312 ELSEIF(KEY3(1:5)=='ktfu8')THEN
313 IKT = 3
314 ELSEIF(KEY3(1:5)=='ktcon')THEN
315 IKT = 4
316 ELSEIF(KEY3(1:5)=='piter')THEN
317 READ(KEY4,'(i5)') IPUPD
318 ELSEIF(KEY3(1:5)=='smdis')THEN
319 ISMDISP = 1
320 ELSEIF(KEY3(1:5)=='solvi')THEN
321 SOLVNFO = 1
322 ELSE
323 READ(KEY3,'(i2)')INSOLV
324 CALL WRIUSC2(IKAD(IKEY)+K,1,KEY0(IKEY))
325 READ(IUSC2,'(a)')TITLE
326 READ(TITLE,*)N_LIM,NITOL,N_TOL
327 IF (NITOL>10) THEN
328 SELECT CASE (NITOL)
329 CASE(12)
330 READ(TITLE,*)N_LIM,NITOL,N_TOLE,N_TOLF
331 CASE(13)
332 READ(TITLE,*)N_LIM,NITOL,N_TOLE,N_TOLU
333 CASE(23)
334 READ(TITLE,*)N_LIM,NITOL,N_TOLF,N_TOLU
335 CASE(123)
336 READ(TITLE,*)N_LIM,NITOL,N_TOLE,N_TOLF,N_TOLU
337 END SELECT
338 ENDIF !(NITOL>10)
339.AND. IF(NITOL==1IRREF==1) IRREF = 0
340 ENDIF
341 ELSEIF(KEY2(1:5)=='sinit')THEN
342 ISIGINI=1
343 ELSEIF(KEY2(1:5)=='lbfgs')THEN
344 READ(KEY3,'(i5)') L_BFGS
345C----------------------------
346C Step Control
347C----------------------------
348 ELSEIF(KEY2=='dtini')THEN
349 CALL WRIUSC2(IKAD(IKEY)+K,1,KEY0(IKEY))
350 READ(IUSC2,*)DT_IMP
351 ELSEIF(KEY2(1:2)=='dt')THEN
352 IF(KEY3(1:4)=='stop')THEN
353 CALL WRIUSC2(IKAD(IKEY)+K,1,KEY0(IKEY))
354 READ(IUSC2,*)DT_MIN,DT_MAX
355C----------------------------
356C-------------fix point for time step----
357C----------------------------
358 ELSEIF(KEY3(1:4)=='fixp')THEN
359 KK =K
360 DO I=1,NBC
361 READ(IUSC1,REC=IKAD(IKEY)+KK,FMT='(a)',ERR=9990)CARTE
362 CALL WRIUSC2(IKAD(IKEY)+KK,1,KEY0(IKEY))
363 NJ = NVAR(CARTE)
364 IF ((NDTFIX+NJ)>100) THEN
365 NJ = 100-NDTFIX
366 WRITE(ISTDO,*)
367 . ' ** warning ** : maximum 100 fix points permitted '
368 ENDIF
369 READ(IUSC2,*,ERR=9990,END=9990)(DTIMPF(NDTFIX+J),J=1,NJ)
370 KK=KK+1
371 NDTFIX = NDTFIX + NJ
372 ENDDO
373 CALL ORDER_DTF(NDTFIX,DTIMPF)
374 ELSE
375 READ(KEY3,'(i2)')IM
376.AND. IF (IDTC>0IM>0) GOTO 9990
377 IDTC=IM
378 IF(IM==1)THEN
379 CALL WRIUSC2(IKAD(IKEY)+K,1,KEY0(IKEY))
380 READ(IUSC2,*)NL_DTP,SCAL_DTP,NL_DTN,SCAL_DTN
381 ELSEIF(IM==2)THEN
382 CALL WRIUSC2(IKAD(IKEY)+K,1,KEY0(IKEY))
383 READ(IUSC2,*)NL_DTP,ALEN0,NL_DTN,SCAL_DTN,SCAL_DTP
384 ELSEIF(IM==3)THEN
385 CALL WRIUSC2(IKAD(IKEY)+K,1,KEY0(IKEY))
386 READ(IUSC2,*)NL_DTP,ALEN0,NL_DTN,SCAL_DTN,SCAL_DTP,IAL_M,
387 . SCAL_RIKS
388 ELSE
389 GOTO 9990
390 ENDIF
391 ENDIF
392C----------------------------
393C NCYCLE stop
394C----------------------------
395 ELSEIF(KEY2=='ncycl')THEN
396 IF(KEY3(1:4)=='stop')THEN
397 CALL WRIUSC2(IKAD(IKEY)+K,1,KEY0(IKEY))
398 READ(IUSC2,*)NCY_MAX
399 ELSE
400 GOTO 9990
401 ENDIF
402C----------------------------
403C interface Control
404C----------------------------
405 ELSEIF(KEY2(1:5)=='inter')THEN
406 IF(KEY3(1:5)=='ttoff')THEN
407 ITTOFF = 1
408 ELSEIF(KEY3(1:5)=='sint7')THEN
409 READ(KEY4,'(i2)')IMP_INT7
410C-----0 nonlinear, 1: linear 2: constant---
411 IMP_INT7= MIN(2,IMP_INT7)
412C---------will be suppressed in the 14.0.210
413 ELSEIF(KEY3(1:5)=='knonl')THEN
414C-----0 nonlinear, 1: linear ----
415 READ(KEY4,'(i2)')IM
416 INTP_C = -IM -1
417C ELSEIF(KEY3(1:5)=='KCOMP'.AND.INTP_C==0)THEN
418 ELSEIF(KEY3(1:5)=='kcomp')THEN
419c INTP_C = 1
420C-----hide option to activate KG in int24 (/IMPLICIT should be defined in Starter)
421 ELSEIF(KEY3(1:4)=='kgon')THEN
422 IIKGOFF = 0
423 ELSE
424 GOTO 9990
425 ENDIF
426C----------------------------
427C R_ref options
428C----------------------------
429 ELSEIF(KEY2(1:4)=='rref')THEN
430 IRREF = 2
431 IF(KEY3(1:3)=='off') THEN
432 IRREF = 0
433 ELSEIF(KEY3(1:5)=='inter')THEN
434C-----0 agressive, 1: moyen ----2--faible 3 non--4 non sauf 1er---
435 READ(KEY4,'(i2)')IM
436 IREFI = IM
437 ELSEIF(KEY3(1:5)=='limit')THEN
438 CALL WRIUSC2(IKAD(IKEY)+K,1,KEY0(IKEY))
439 READ(IUSC2,*)RF_MIN,RF_MAX
440 ENDIF
441C----------------------------
442C divergence criteria
443C----------------------------
444 ELSEIF(KEY2(1:5)=='diver')THEN
445 IF(KEY3(1:3)=='tol')THEN
446 CALL WRIUSC2(IKAD(IKEY)+K,1,KEY0(IKEY))
447 READ(IUSC2,*)TOL_DIV
448C-----num. of diver---
449 ELSE
450 READ(KEY3,'(i2)')IM
451 NDIVER = IM
452 IF (NDIVER ==0) NDIVER=-1
453 END IF
454C----------------------------
455C Geometrical stifness
456C----------------------------
457 ELSEIF(KEY2(1:5)=='gstif')THEN
458 IF(KEY3(1:3)=='off')IKG=0
459C----------------------------
460C Geometrical stifness
461C----------------------------
462 ELSEIF(KEY2(1:5)=='pstif')THEN
463 IF(KEY3(1:3)=='off') IKPRES=0
464C----------------------------
465C buckling analysis
466C----------------------------
467 ELSEIF(KEY2=='buckl')THEN
468 READ(KEY3,'(i2)')IBUCKL
469 IF (IBUCKL==0) THEN
470 WRITE(ISTDO,*) ' ** error ** : keyword /impl/buckl obsolete ',
471 . 'using /impl/buckl/1 or /impl/buckl/2'
472 GOTO 9990
473 ENDIF
474 IBUCKL = IBUCKL-1
475 CALL WRIUSC2(IKAD(IKEY)+K,1,KEY0(IKEY))
476 READ(IUSC2,*) EMIN_B, EMAX_B, NBUCK, MSGL_B, MAXSET_B, SHIFT_B
477 IF (SHIFT_B==ZERO) SHIFT_B=EM02
478 SHFTBUCK = SHIFT_B
479 IF (MAXSET_B==0) MAXSET_B=8
480 BNITER=300
481 BINCV=4
482 BMAXNCV=16
483c BMAXNCV=MAX(BINCV,BMAXNCV)
484 BIPRI =MSGL_B
485 BISOLV=1
486C
487 ELSEIF(KEY2(1:5)=='autos')THEN
488 IF(KEY3(1:3)=='off')THEN
489 IAUTSPC=0
490 ELSEIF(KEY3(1:3)=='all')THEN
491 IAUTSPC=2
492 ENDIF
493C----------------------------
494C line_search option
495C--------0=3, 1:energy 2: force --3:AUTO (old)------------------
496 ELSEIF(KEY2(1:5)=='lsear')THEN
497 IF(KEY3(1:3)=='off')THEN
498 ILINE_S = 100
499 ELSE
500 READ(KEY3,'(i2)')ILINE_S
501 CALL WRIUSC2(IKAD(IKEY)+K,1,KEY0(IKEY))
502 READ(IUSC2,*)NLS_LIM,LS_TOL
503 ENDIF
504C----------------------------
505C projection for warped shell elements
506C--------0=no proj but keep explicit part 1:doing -1 :no proj, neither for explicit---
507 ELSEIF(KEY2(1:5)=='shpof')THEN
508 IKPROJ=-1
509C-----------become default-after-----
510 ELSEIF(KEY2(1:5)=='shpon')THEN
511 IKPROJ=1
512C----------------------------
513C OLD CONTROL OPTIONS
514C----------------------------
515 ELSEIF(KEY2(1:5)=='contr')THEN
516 IF(KEY3(1:2)=='dt')THEN
517 IF(KEY4(1:4)=='stop')THEN
518 CALL WRIUSC2(IKAD(IKEY)+K,1,KEY0(IKEY))
519 READ(IUSC2,*)DT_MIN,DT_MAX
520 ELSE
521 READ(KEY4,'(i2)')IM
522 IDTC=IM
523 IF(IM==1)THEN
524 CALL WRIUSC2(IKAD(IKEY)+K,1,KEY0(IKEY))
525 READ(IUSC2,*)NL_DTP,SCAL_DTP,NL_DTN,SCAL_DTN
526 ELSEIF(IM==2)THEN
527 CALL WRIUSC2(IKAD(IKEY)+K,1,KEY0(IKEY))
528 READ(IUSC2,*)NL_DTP,ALEN0,NL_DTN,SCAL_DTN,SCAL_DTP
529 ENDIF
530 ENDIF
531 ELSEIF(KEY3(1:4)=='shel')THEN
532C----------------------------
533C Fictif stifness of Mzz for shell
534C----------------------------
535 CALL WRIUSC2(IKAD(IKEY)+K,1,KEY0(IKEY))
536 READ(IUSC2,*)KZ_TOL
537 ELSEIF(KEY3(1:5)=='inter')THEN
538C----------------------------
539C stifness factor for interface
540C----------------------------
541 CALL WRIUSC2(IKAD(IKEY)+K,1,KEY0(IKEY))
542 READ(IUSC2,*)SK_INT
543 ENDIF
544C----------------------------
545C hide options
546C----------------------------
547 ELSEIF(KEY2(1:5)=='prtol')THEN
548 CALL WRIUSC2(IKAD(IKEY)+K,1,KEY0(IKEY))
549 READ(IUSC2,*)D_TOL
550 ELSEIF(KEY2(1:4)=='nexp')THEN
551 READ(KEY3,'(i5)')NEXP
552 ELSEIF(KEY2=='debug')THEN
553 impdeb=1
554 CALL wriusc2(ikad(ikey)+k,1,key0(ikey))
555 READ(iusc2,*)ndeb0,ndeb1
556 IF(ndeb0/=0)ndeb0 = ndeb0 + 1
557 ndeb1=max(ndeb0,ndeb1+1)
558 ELSEIF(key2(1:3)=='DEL')THEN
559 IF(key3(1:5)=='RBODY')THEN
560 imp_rby=1
561 ELSEIF(key3(1:5)=='INTER')THEN
562 imp_int=1
563 ENDIF
564 ELSEIF(key2(1:5)=='ITRBY')THEN
565C-------max iter for secnd dis calculation with finite rotation---
566 READ(key3,'(I3)')itrmax
567 ELSEIF(key2(1:4)=='LRIG')THEN
568 imp_lr = 1
569 ELSE
570 GOTO 9990
571 ENDIF
572 k=k+nbc
573 IF(ikad(ikey)+k/=ikad(ikey+1))GO TO 1160
574 IF (iparit/=0) THEN
575 iparit=0
576 ikg=ikg+5
577 ENDIF
578 ENDIF
579C
580 RETURN
581C
582 9990 CONTINUE
583 CALL ancmsg(msgid=73,anmode=aninfo,
584 . c1=key0(ikey))
585 CALL arret(0)
#define max(a, b)
Definition macros.h:21
integer b_mcore
integer msg_lvl
integer b_order
integer, dimension(4) e_ref
integer, parameter ncharline100
subroutine ancmsg(msgid, msgtype, anmode, i1, i2, i3, i4, i5, i6, i7, i8, i9, i10, i11, i12, i13, i14, i15, i16, i17, i18, i19, i20, r1, r2, r3, r4, r5, r6, r7, r8, r9, c1, c2, c3, c4, c5, c6, c7, c8, c9, prmode)
Definition message.F:889
subroutine arret(nn)
Definition arret.F:87
subroutine wriusc2(irec, nbc, key0)
Definition wriusc2.F:60

◆ order_dtf()

subroutine order_dtf ( integer n,
rc )

Definition at line 592 of file freimpl.F.

593C----6---------------------------------------------------------------7---------8
594C I m p l i c i t T y p e s
595C-----------------------------------------------
596#include "implicit_f.inc"
597C-----------------------------------------------------------------
598C D u m m y A r g u m e n t s
599C-----------------------------------------------
600 INTEGER N
601 my_real
602 . rc(*)
603C-----------------------------------------------
604C L o c a l V a r i a b l e s
605C-----------------------------------------------
606 INTEGER I,J,II,NN
607 my_real
608 . s(n),smin
609C
610 IF (n==0) RETURN
611C-----en ordre croisante-----
612 nn =0
613 DO i =1,n
614 IF (rc(i)>zero) THEN
615 nn = nn +1
616 s(nn)= rc(i)
617 ENDIF
618 ENDDO
619 n= nn
620 DO i =1,n
621 smin=s(i)
622 ii=i
623 DO j =i+1,n
624 IF (s(j)<smin) THEN
625 ii=j
626 smin = s(j)
627 ENDIF
628 ENDDO
629 IF (ii/=i) THEN
630 smin =s(i)
631 s(i)=s(ii)
632 s(ii)=smin
633 ENDIF
634 rc(i) = s(i)
635 ENDDO
636C----6---------------------------------------------------------------7---------8
637 RETURN
#define my_real
Definition cppsort.cpp:32