OpenRadioss 2025.1.11
OpenRadioss project
Loading...
Searching...
No Matches
nl_solv.F File Reference
#include "implicit_f.inc"
#include "com01_c.inc"
#include "com04_c.inc"
#include "com08_c.inc"
#include "impl1_c.inc"
#include "impl2_c.inc"
#include "units_c.inc"
#include "task_c.inc"

Go to the source code of this file.

Functions/Subroutines

subroutine nl_solv (nddl, iddl, ndof, ikc, d, dr, nnz, iadk, jdik, diag_k, lt_k, f, nddli, iadi, jdii, diag_i, lt_i, itok, iadm, jdim, diag_m, lt_m, r02, dd, ddr, itask0, it, itc, ru0, rold, idiv, inprint, icprec, istop, e02, de0, eimp, inloc, nddl0, ls, u02, gap, itab, fr_elem, iad_elem, w_ddl, a, ar, v, ms, x, ipari, intbuf_tab, num_imp, ns_imp, ne_imp, nsrem, nsl, icont, graphe, fac_k, ipiv_k, nk, nmonv, imonv, monvol, igrsurf, fr_mv, volmon, ibfv, skew, xframe, mumps_par, cddlp, ind_imp, nbintc, intlist, newfront, isendto, irecvfrom, irbe3, lrbe3, ndiv, icont0, isign, fext, dg, dgr, dg0, dgr0, rfext, ls1, nodft, nodlt, irbe2, lrbe2, idiv0, relres, anew_stif)
subroutine crit_ite (it, ur, rr, er, ndiv, tol)
subroutine line_s (r0, s0, r1, s, dtol, idiv, iint, prec)
subroutine get_max (n, d, vmax, i)
subroutine line_s1 (e1, s, ep, sp, en, sn, idiv, dtol, icont, icont0, iriks)
subroutine al_constraint1 (nddl0, nddl, iddl, ndof, ikc, dd, ddr, dg, dgr, di, dir, w_ddl, l_a, lamda, sw2, ier)
subroutine al_constraint2 (nddl0, nddl, iddl, ndof, ikc, dd, ddr, dg, dgr, di, dir, w_ddl, l_a, lamda, sw2)
subroutine b_correct (nddl, f, fext, dlamda)
subroutine b_correct_h (f_ddl, l_ddl, f, fext, dlamda)
subroutine al_constrainth1 (nddl0, nddl, iddl, ndof, ikc, dd, ddr, dg, dgr, di, dir, w_ddl, l_a, lamda, sw2, ft_n, lt_n, f_ddl, l_ddl, itask, ier)
subroutine al_constrainth2 (nddl0, nddl, iddl, ndof, ikc, dd, ddr, dg, dgr, di, dir, w_ddl, l_a, lamda, sw2, ft_n, lt_n, f_ddl, l_ddl, itask)
subroutine b_correct_hp (nddl, f, fext, dlamda)
subroutine al_constraint1_hp (nddl0, nddl, iddl, ndof, ikc, dd, ddr, dg, dgr, di, dir, w_ddl, l_a, lamda, sw2, ier)
subroutine al_constraint2_hp (nddl0, nddl, iddl, ndof, ikc, dd, ddr, dg, dgr, di, dir, w_ddl, l_a, lamda, sw2)

Function/Subroutine Documentation

◆ al_constraint1()

subroutine al_constraint1 ( integer nddl0,
integer nddl,
integer, dimension(*) iddl,
integer, dimension(*) ndof,
integer, dimension(*) ikc,
dd,
ddr,
dg,
dgr,
di,
dir,
integer, dimension(*) w_ddl,
l_a,
lamda,
sw2,
integer ier )

Definition at line 1115 of file nl_solv.F.

1119C-----------------------------------------------
1120C I m p l i c i t T y p e s
1121C-----------------------------------------------
1122#include "implicit_f.inc"
1123C-----------------------------------------------
1124C C o m m o n B l o c k s
1125C-----------------------------------------------
1126#include "com04_c.inc"
1127C-----------------------------------------------
1128C D u m m y A r g u m e n t s
1129C-----------------------------------------------
1130 INTEGER NDDL,NDDL0,IDDL(*) ,NDOF(*) ,IKC(*),W_DDL(*) ,IER
1131C REAL
1132 my_real
1133 . lamda,dg(3,*),dgr(3,*),di(3,*),dir(3,*),dd(3,*),ddr(3,*),
1134 . l_a,sw2
1135C-----------------------------------------------
1136C L o c a l V a r i a b l e s
1137C-----------------------------------------------
1138 INTEGER I
1139 my_real
1140 . a,b,c,ui2,uiud,uiut,ut2,lam1,lam2,s,th1,th2,fac
1141C-----------------------------------------------
1142 ier=0
1143 CALL produt_u(nddl0 ,nddl ,iddl ,ndof ,ikc ,
1144 . dg ,dgr ,ut2 ,w_ddl )
1145 a = ut2+sw2
1146C-------DD< DI+DD
1147 CALL integrator2(1,numnod,iddl,ndof,ikc,dd,ddr,di,dir)
1148 CALL produt_u2(nddl0 ,nddl ,iddl ,ndof ,ikc ,
1149 . dg ,dgr ,dd ,ddr ,uiut ,
1150 . w_ddl )
1151 b= two*(uiut+sw2*lamda)
1152 CALL produt_u(nddl0 ,nddl ,iddl ,ndof ,ikc ,
1153 . dd ,ddr ,ui2 ,w_ddl )
1154 c = ui2-l_a*l_a+sw2*lamda
1155 s= b*b-four*a*c
1156 IF (s>=0) THEN
1157 s=sqrt(s)
1158 lam1=half*(-b+s)/a
1159 lam2=half*(-b-s)/a
1160 CALL produt_u2(nddl0 ,nddl ,iddl ,ndof ,ikc ,
1161 . dd ,ddr ,di ,dir ,uiud ,
1162 . w_ddl )
1163 CALL produt_u2(nddl0 ,nddl ,iddl ,ndof ,ikc ,
1164 . dg ,dgr ,di ,dir ,uiut ,
1165 . w_ddl )
1166 th1= uiud+lam1*uiut
1167 th2= uiud+lam2*uiut
1168 IF (th1>zero.AND.th2<zero) THEN
1169 lamda =lam1
1170 ELSEIF (th1<zero.AND.th2>zero) THEN
1171 lamda =lam2
1172 ELSEIF (th1>zero.AND.th2>zero) THEN
1173 lamda =min(lam1,lam2)
1174 ELSE
1175 ier=1
1176 ENDIF
1177 ELSE
1178 ier=1
1179 ENDIF
1180 fac =-one
1181 CALL frac_dd(1,numnod,iddl,ndof,ikc,dd,ddr,di,dir,fac)
1182C
1183 RETURN
my_real
#define my_real
Definition cppsort.cpp:32
integrator2
subroutine integrator2(nodft, nodlt, iddl, ndof, ikc, d, dr, dd, ddr)
Definition integrator.F:207
frac_dd
subroutine frac_dd(nodft, nodlt, iddl, ndof, ikc, d, dr, dd, ddr, fac)
Definition integrator.F:295
min
#define min(a, b)
Definition macros.h:20
produt_u
subroutine produt_u(nddl0, nddl, iddl, ndof, ikc, dd, ddr, norm2, w_imp)
Definition produt_v.F:224
produt_u2
subroutine produt_u2(nddl0, nddl, iddl, ndof, ikc, d1, d1r, d2, d2r, norm2, w_imp)
Definition produt_v.F:262

◆ al_constraint1_hp()

subroutine al_constraint1_hp ( integer nddl0,
integer nddl,
integer, dimension(*) iddl,
integer, dimension(*) ndof,
integer, dimension(*) ikc,
dd,
ddr,
dg,
dgr,
di,
dir,
integer, dimension(*) w_ddl,
l_a,
lamda,
sw2,
integer ier )

Definition at line 1459 of file nl_solv.F.

1463C-----------------------------------------------
1464C I m p l i c i t T y p e s
1465C-----------------------------------------------
1466#include "implicit_f.inc"
1467C-----------------------------------------------
1468C D u m m y A r g u m e n t s
1469C-----------------------------------------------
1470 INTEGER NDDL,NDDL0,IDDL(*) ,NDOF(*) ,IKC(*),W_DDL(*) ,IER,
1471 . FT_N ,LT_N
1472C REAL
1473 my_real
1474 . lamda,dg(3,*),dgr(3,*),di(3,*),dir(3,*),dd(3,*),ddr(3,*),
1475 . l_a,sw2
1476C-----------------------------------------------
1477C L o c a l V a r i a b l e s
1478C-----------------------------------------------
1479 INTEGER I
1480 my_real
1481 . a,b,c,ui2,uiud,uiut,ut2,lam1,lam2,s,th1,th2,fac
1482C-----------------------------------------------
1483 ier=0
1484 CALL produt_uhp(nddl0 ,nddl ,iddl ,ndof ,ikc ,
1485 . dg ,dgr ,ut2 ,w_ddl )
1486C---------------------
1487 a = ut2+sw2
1488C-------DD< DI+DD
1489 CALL integrator2_hp(iddl,ndof,ikc,dd,ddr,di,dir)
1490C---------------------
1491 CALL produt_uhp2(nddl0 ,nddl ,iddl ,ndof ,ikc ,
1492 . dg ,dgr ,dd ,ddr ,uiut ,
1493 . w_ddl )
1494 CALL produt_uhp(nddl0 ,nddl ,iddl ,ndof ,ikc ,
1495 . dd ,ddr ,ui2 ,w_ddl )
1496C---------------------
1497 b= two*(uiut+sw2*lamda)
1498 c = ui2-l_a*l_a+sw2*lamda
1499 s= b*b-four*a*c
1500 IF (s>=0) THEN
1501 s=sqrt(s)
1502 lam1=half*(-b+s)/a
1503 lam2=half*(-b-s)/a
1504 CALL produt_uhp2(nddl0 ,nddl ,iddl ,ndof ,ikc ,
1505 . dd ,ddr ,di ,dir ,uiud ,
1506 . w_ddl )
1507 CALL produt_uhp2(nddl0 ,nddl ,iddl ,ndof ,ikc ,
1508 . dg ,dgr ,di ,dir ,uiut ,
1509 . w_ddl )
1510C---------------------
1511 th1= uiud+lam1*uiut
1512 th2= uiud+lam2*uiut
1513 IF (th1>zero.AND.th2<zero) THEN
1514 lamda =lam1
1515 ELSEIF (th1<zero.AND.th2>zero) THEN
1516 lamda =lam2
1517 ELSEIF (th1>zero.AND.th2>zero) THEN
1518 lamda =min(lam1,lam2)
1519 ELSE
1520 ier=1
1521 ENDIF
1522 ELSE
1523 ier=1
1524 ENDIF
1525 fac =-one
1526 CALL frac_dd_hp(iddl,ndof,ikc,dd,ddr,di,dir,fac)
1527C
1528 RETURN
frac_dd_hp
subroutine frac_dd_hp(iddl, ndof, ikc, d, dr, dd, ddr, fac)
Definition integrator.F:618
integrator2_hp
subroutine integrator2_hp(iddl, ndof, ikc, d, dr, dd, ddr)
Definition integrator.F:513
produt_uhp
subroutine produt_uhp(nddl0, nddl, iddl, ndof, ikc, dd, ddr, norm2, w_imp)
Definition produt_v.F:3359
produt_uhp2
subroutine produt_uhp2(nddl0, nddl, iddl, ndof, ikc, d1, d1r, d2, d2r, norm2, w_imp)
Definition produt_v.F:3398

◆ al_constraint2()

subroutine al_constraint2 ( integer nddl0,
integer nddl,
integer, dimension(*) iddl,
integer, dimension(*) ndof,
integer, dimension(*) ikc,
dd,
ddr,
dg,
dgr,
di,
dir,
integer, dimension(*) w_ddl,
l_a,
lamda,
sw2 )

Definition at line 1191 of file nl_solv.F.

1194C-----------------------------------------------
1195C I m p l i c i t T y p e s
1196C-----------------------------------------------
1197#include "implicit_f.inc"
1198C-----------------------------------------------
1199C D u m m y A r g u m e n t s
1200C-----------------------------------------------
1201 INTEGER NDDL,NDDL0,IDDL(*) ,NDOF(*) ,IKC(*),W_DDL(*)
1202 my_real
1203 . lamda,dg(3,*),dgr(3,*),di(3,*),dir(3,*),dd(3,*),ddr(3,*),
1204 . l_a,sw2
1205C-----------------------------------------------
1206C L o c a l V a r i a b l e s
1207C-----------------------------------------------
1208 INTEGER I
1209 my_real
1210 . gi,ui2,uiud,uiut,s
1211C-----------------------------------------------
1212 CALL produt_u(nddl0 ,nddl ,iddl ,ndof ,ikc ,
1213 . di ,dir ,ui2 ,w_ddl )
1214 gi = sqrt(ui2)
1215 CALL produt_u2(nddl0 ,nddl ,iddl ,ndof ,ikc ,
1216 . dg ,dgr ,di ,dir ,uiut ,
1217 . w_ddl )
1218 CALL produt_u2(nddl0 ,nddl ,iddl ,ndof ,ikc ,
1219 . dd ,ddr ,di ,dir ,uiud ,
1220 . w_ddl )
1221 s = sw2*lamda
1222 lamda = ((l_a-gi)*gi-uiud)/(uiut+s)
1223C
1224 RETURN

◆ al_constraint2_hp()

subroutine al_constraint2_hp ( integer nddl0,
integer nddl,
integer, dimension(*) iddl,
integer, dimension(*) ndof,
integer, dimension(*) ikc,
dd,
ddr,
dg,
dgr,
di,
dir,
integer, dimension(*) w_ddl,
l_a,
lamda,
sw2 )

Definition at line 1538 of file nl_solv.F.

1541C-----------------------------------------------
1542C I m p l i c i t T y p e s
1543C-----------------------------------------------
1544#include "implicit_f.inc"
1545C-----------------------------------------------
1546C D u m m y A r g u m e n t s
1547C-----------------------------------------------
1548 INTEGER NDDL,NDDL0,IDDL(*) ,NDOF(*) ,IKC(*),W_DDL(*)
1549C REAL
1550 my_real
1551 . lamda,dg(3,*),dgr(3,*),di(3,*),dir(3,*),dd(3,*),ddr(3,*),
1552 . l_a,sw2
1553C-----------------------------------------------
1554C L o c a l V a r i a b l e s
1555C-----------------------------------------------
1556 INTEGER I
1557 my_real
1558 . gi,ui2,uiud,uiut,s
1559C-----------------------------------------------
1560 CALL produt_uhp(nddl0 ,nddl ,iddl ,ndof ,ikc ,
1561 . di ,dir ,ui2 ,w_ddl )
1562 CALL produt_uhp2(nddl0 ,nddl ,iddl ,ndof ,ikc ,
1563 . dg ,dgr ,di ,dir ,uiut ,
1564 . w_ddl )
1565 CALL produt_uhp2(nddl0 ,nddl ,iddl ,ndof ,ikc ,
1566 . dd ,ddr ,di ,dir ,uiud ,
1567 . w_ddl )
1568C---------------------
1569 gi = sqrt(ui2)
1570 s = sw2*lamda
1571 lamda = ((l_a-gi)*gi-uiud)/(uiut+s)
1572C
1573 RETURN

◆ al_constrainth1()

subroutine al_constrainth1 ( integer nddl0,
integer nddl,
integer, dimension(*) iddl,
integer, dimension(*) ndof,
integer, dimension(*) ikc,
dd,
ddr,
dg,
dgr,
di,
dir,
integer, dimension(*) w_ddl,
l_a,
lamda,
sw2,
integer ft_n,
integer lt_n,
integer f_ddl,
integer l_ddl,
integer itask,
integer ier )

Definition at line 1286 of file nl_solv.F.

1291C-----------------------------------------------
1292C I m p l i c i t T y p e s
1293C-----------------------------------------------
1294#include "implicit_f.inc"
1295C-----------------------------------------------
1296C D u m m y A r g u m e n t s
1297C-----------------------------------------------
1298 INTEGER NDDL,NDDL0,IDDL(*) ,NDOF(*) ,IKC(*),W_DDL(*) ,IER,
1299 . FT_N ,LT_N ,F_DDL ,L_DDL ,ITASK
1300C REAL
1301 my_real
1302 . lamda,dg(3,*),dgr(3,*),di(3,*),dir(3,*),dd(3,*),ddr(3,*),
1303 . l_a,sw2
1304C-----------------------------------------------
1305C L o c a l V a r i a b l e s
1306C-----------------------------------------------
1307 INTEGER I
1308 my_real
1309 . a,b,c,ui2,uiud,uiut,ut2,lam1,lam2,s,th1,th2,fac
1310C-----------------------------------------------
1311 ier=0
1312 CALL produt_uh(nddl0 ,nddl ,iddl ,ndof ,ikc ,
1313 . dg ,dgr ,ut2 ,w_ddl ,f_ddl ,
1314 . l_ddl ,itask )
1315C----------------------
1316 CALL my_barrier
1317C---------------------
1318 a = ut2+sw2
1319C-------DD< DI+DD
1320 CALL integrator2(ft_n,lt_n,iddl,ndof,ikc,dd,ddr,di,dir)
1321C----------------------
1322 CALL my_barrier
1323C---------------------
1324 CALL produt_uh2(nddl0 ,nddl ,iddl ,ndof ,ikc ,
1325 . dg ,dgr ,dd ,ddr ,uiut ,
1326 . w_ddl ,f_ddl ,l_ddl ,itask )
1327 CALL produt_uh(nddl0 ,nddl ,iddl ,ndof ,ikc ,
1328 . dd ,ddr ,ui2 ,w_ddl ,f_ddl ,
1329 . l_ddl ,itask )
1330C----------------------
1331 CALL my_barrier
1332C---------------------
1333 b= two*(uiut+sw2*lamda)
1334 c = ui2-l_a*l_a+sw2*lamda
1335 s= b*b-four*a*c
1336 IF (s>=0) THEN
1337 s=sqrt(s)
1338 lam1=half*(-b+s)/a
1339 lam2=half*(-b-s)/a
1340 CALL produt_uh2(nddl0 ,nddl ,iddl ,ndof ,ikc ,
1341 . dd ,ddr ,di ,dir ,uiud ,
1342 . w_ddl ,f_ddl ,l_ddl ,itask )
1343 CALL produt_uh2(nddl0 ,nddl ,iddl ,ndof ,ikc ,
1344 . dg ,dgr ,di ,dir ,uiut ,
1345 . w_ddl ,f_ddl ,l_ddl ,itask )
1346C----------------------
1347 CALL my_barrier
1348C---------------------
1349 th1= uiud+lam1*uiut
1350 th2= uiud+lam2*uiut
1351 IF (th1>zero.AND.th2<zero) THEN
1352 lamda =lam1
1353 ELSEIF (th1<zero.AND.th2>zero) THEN
1354 lamda =lam2
1355 ELSEIF (th1>zero.AND.th2>zero) THEN
1356 lamda =min(lam1,lam2)
1357 ELSE
1358 ier=1
1359 ENDIF
1360 ELSE
1361 ier=1
1362 ENDIF
1363 fac =-one
1364 CALL frac_dd(ft_n,lt_n,iddl,ndof,ikc,dd,ddr,di,dir,fac)
1365C
1366 RETURN
produt_uh2
subroutine produt_uh2(nddl0, nddl, iddl, ndof, ikc, d1, d1r, d2, d2r, norm2, w_imp, f_ddl, l_ddl, itask)
Definition produt_v.F:1741
produt_uh
subroutine produt_uh(nddl0, nddl, iddl, ndof, ikc, dd, ddr, norm2, w_imp, f_ddl, l_ddl, itask)
Definition produt_v.F:1687
my_barrier
subroutine my_barrier
Definition machine.F:31

◆ al_constrainth2()

subroutine al_constrainth2 ( integer nddl0,
integer nddl,
integer, dimension(*) iddl,
integer, dimension(*) ndof,
integer, dimension(*) ikc,
dd,
ddr,
dg,
dgr,
di,
dir,
integer, dimension(*) w_ddl,
l_a,
lamda,
sw2,
integer ft_n,
integer lt_n,
integer f_ddl,
integer l_ddl,
integer itask )

Definition at line 1375 of file nl_solv.F.

1379C-----------------------------------------------
1380C I m p l i c i t T y p e s
1381C-----------------------------------------------
1382#include "implicit_f.inc"
1383C-----------------------------------------------
1384C D u m m y A r g u m e n t s
1385C-----------------------------------------------
1386 INTEGER NDDL,NDDL0,IDDL(*) ,NDOF(*) ,IKC(*),W_DDL(*),
1387 . FT_N ,LT_N ,F_DDL ,L_DDL ,ITASK
1388C REAL
1389 my_real
1390 . lamda,dg(3,*),dgr(3,*),di(3,*),dir(3,*),dd(3,*),ddr(3,*),
1391 . l_a,sw2
1392C-----------------------------------------------
1393C L o c a l V a r i a b l e s
1394C-----------------------------------------------
1395 INTEGER I
1396 my_real
1397 . gi,ui2,uiud,uiut,s
1398C-----------------------------------------------
1399 CALL produt_uh(nddl0 ,nddl ,iddl ,ndof ,ikc ,
1400 . di ,dir ,ui2 ,w_ddl ,f_ddl ,
1401 . l_ddl ,itask )
1402 CALL produt_uh2(nddl0 ,nddl ,iddl ,ndof ,ikc ,
1403 . dg ,dgr ,di ,dir ,uiut ,
1404 . w_ddl ,f_ddl ,l_ddl ,itask )
1405 CALL produt_uh2(nddl0 ,nddl ,iddl ,ndof ,ikc ,
1406 . dd ,ddr ,di ,dir ,uiud ,
1407 . w_ddl ,f_ddl ,l_ddl ,itask )
1408C----------------------
1409 CALL my_barrier
1410C---------------------
1411 gi = sqrt(ui2)
1412 s = sw2*lamda
1413 lamda = ((l_a-gi)*gi-uiud)/(uiut+s)
1414C
1415 RETURN

◆ b_correct()

subroutine b_correct ( integer nddl,
f,
fext,
dlamda )

Definition at line 1229 of file nl_solv.F.

1230C-----------------------------------------------
1231C I m p l i c i t T y p e s
1232C-----------------------------------------------
1233#include "implicit_f.inc"
1234C-----------------------------------------------
1235C D u m m y A r g u m e n t s
1236C-----------------------------------------------
1237 INTEGER NDDL
1238 my_real
1239 . dlamda,f(*),fext(*)
1240C-----------------------------------------------
1241C L o c a l V a r i a b l e s
1242C-----------------------------------------------
1243 INTEGER I
1244C-----------------------------------------------
1245 DO i =1,nddl
1246 f(i)=f(i)+dlamda*fext(i)
1247 ENDDO
1248C
1249 RETURN

◆ b_correct_h()

subroutine b_correct_h ( integer f_ddl,
integer l_ddl,
f,
fext,
dlamda )

Definition at line 1254 of file nl_solv.F.

1255C-----------------------------------------------
1256C I m p l i c i t T y p e s
1257C-----------------------------------------------
1258#include "implicit_f.inc"
1259C-----------------------------------------------
1260C D u m m y A r g u m e n t s
1261C-----------------------------------------------
1262 INTEGER F_DDL ,L_DDL
1263C REAL
1264 my_real
1265 . dlamda,f(*),fext(*)
1266C-----------------------------------------------
1267C L o c a l V a r i a b l e s
1268C-----------------------------------------------
1269 INTEGER I
1270C-----------------------------------------------
1271 DO i =f_ddl ,l_ddl
1272 f(i)=f(i)+dlamda*fext(i)
1273 ENDDO
1274C
1275 RETURN

◆ b_correct_hp()

subroutine b_correct_hp ( integer nddl,
f,
fext,
dlamda )

Definition at line 1422 of file nl_solv.F.

1423C-----------------------------------------------
1424C I m p l i c i t T y p e s
1425C-----------------------------------------------
1426#include "implicit_f.inc"
1427C-----------------------------------------------
1428C D u m m y A r g u m e n t s
1429C-----------------------------------------------
1430 INTEGER NDDL
1431C REAL
1432 my_real
1433 . dlamda,f(*),fext(*)
1434C-----------------------------------------------
1435C L o c a l V a r i a b l e s
1436C-----------------------------------------------
1437 INTEGER F_DDL ,L_DDL ,ITSK
1438 INTEGER I
1439C-----------------------------------------------
1440!$OMP PARALLEL PRIVATE(ITSK,F_DDL ,L_DDL,I)
1441 CALL imp_smpini(itsk ,f_ddl ,l_ddl ,nddl )
1442 DO i =f_ddl ,l_ddl
1443 f(i)=f(i)+dlamda*fext(i)
1444 ENDDO
1445!$OMP END PARALLEL
1446C
1447 RETURN
imp_smpini
subroutine imp_smpini(itsk, n1ftsk, n1ltsk, n1)
Definition imp_solv.F:6895

◆ crit_ite()

subroutine crit_ite ( integer it,
ur,
rr,
er,
integer ndiv,
tol )

Definition at line 808 of file nl_solv.F.

809C-----------------------------------------------
810C I m p l i c i t T y p e s
811C-----------------------------------------------
812#include "implicit_f.inc"
813C-----------------------------------------------
814C C o m m o n B l o c k s
815C-----------------------------------------------
816#include "impl1_c.inc"
817#include "impl2_c.inc"
818C-----------------------------------------------
819C D u m m y A r g u m e n t s
820C-----------------------------------------------
821 INTEGER IT,NDIV
822C REAL
823 my_real
824 . ur,rr,er,tol
825C-----------------------------------------------
826C L o c a l V a r i a b l e s
827C-----------------------------------------------
828 my_real
829 . err,tolr,told
830C--------cas particulier---------------------------------------
831C TOLD=EP04
832C IF (ILINE_S==1) TOLD=EP03
833C IF (ISMDISP==1) TOLD=EP10
834 told = tol_div
835 tolr = one+tol
836C--------1:E,2:r,3:E+U,4:r+U,5,r+E,6,r+E+U-----
837C--new: 1:E,2:r,3:U,12:E+r,13,E+U,23,r+U,123,r+E+U
838 IF (nitol==1) THEN
839 err=abs(er)/tol
840 IF (it==1.OR.rr>=one) err=err+tolr
841 ELSEIF (nitol==2) THEN
842 err=rr/tol
843 ELSEIF (nitol==3) THEN
844 err=ur/tol
845C ERR=MAX(SQRT(ABS(ER)),UR)
846 IF (rr>=one) err=err+tolr
847 ELSEIF (nitol==12) THEN
848 err=max(rr/n_tolf,abs(er)/n_tole)
849c ERR=MAX(RR,UR)
850 ELSEIF (nitol==13) THEN
851 err=max(ur/n_tolu,abs(er)/n_tole)
852c ERR=MAX(SQRT(ABS(ER)),RR)
853 ELSEIF (nitol==23) THEN
854 err=max(ur/n_tolu,rr/n_tolf)
855c ERR=MAX(SQRT(ABS(ER)),RR,UR)
856 ELSEIF (nitol==123) THEN
857 err=max(ur/n_tolu,rr/n_tolf,abs(er)/n_tole)
858 ENDIF
859C
860 IF (err<=one) THEN
861C IF (ERR<=TOL) THEN
862 imconv=1
863!sb ELSEIF (UR<MIN(EM05,TOL).AND.RR<EM01) THEN
864 ELSEIF ((n_lim /= 1.OR.isprb /= 0).AND.ur<min(em05,tol).AND.rr<em01) THEN
865 imconv=1
866 ELSEIF (imconv<0) THEN
867C--------line-search-----
868 ELSE
869 imconv=0
870C------divergence-----------
871!sb
872 IF(n_lim /= 1 .OR. isprb /= 0) THEN
873 IF (rr>tolr) THEN
874 IF (it==1) THEN
875 imconv=-3
876 IF (ndiv<ndiver) imconv=0
877 ELSE
878 imconv=-2
879 IF ((ndiv<ndiver.OR.ilast==1).AND.rr<told) imconv=0
880 ENDIF
881 ENDIF
882 ENDIF
883
884 IF (idtc>0.AND.it>nl_dtn) imconv=-2
885 ENDIF
886C-----------------------------------------------
887 RETURN
max
#define max(a, b)
Definition macros.h:21

◆ get_max()

subroutine get_max ( integer n,
d,
vmax,
integer i )

Definition at line 961 of file nl_solv.F.

962C-----------------------------------------------
963C I m p l i c i t T y p e s
964C-----------------------------------------------
965#include "implicit_f.inc"
966C-----------------------------------------------
967C C o m m o n B l o c k s
968C-----------------------------------------------
969#include "com01_c.inc"
970C-----------------------------------------------
971C D u m m y A r g u m e n t s
972C-----------------------------------------------
973 INTEGER N,I
974 my_real
975 . d(3,*),vmax
976C-----------------------------------------------
977C L o c a l V a r i a b l e s
978C-----------------------------------------------
979 INTEGER K
980 my_real
981 . t
982C-----------------------------------------------
983 vmax=zero
984 i=1
985 DO k=1,n
986 t=max(abs(d(1,k)),abs(d(2,k)),abs(d(3,k)))
987 IF (t>vmax) THEN
988 vmax = t
989 i=k
990 ENDIF
991 ENDDO
992 IF (nspmd>1) CALL spmd_max_s(vmax)
993C
994 RETURN
spmd_max_s
subroutine spmd_max_s(s)
Definition imp_spmd.F:1195

◆ line_s()

subroutine line_s ( r0,
s0,
r1,
s,
dtol,
integer idiv,
integer iint,
prec )

Definition at line 894 of file nl_solv.F.

895C-----------------------------------------------
896C I m p l i c i t T y p e s
897C-----------------------------------------------
898#include "implicit_f.inc"
899C-----------------------------------------------
900C C o m m o n B l o c k s
901C-----------------------------------------------
902#include "impl1_c.inc"
903C-----------------------------------------------
904C D u m m y A r g u m e n t s
905C-----------------------------------------------
906 INTEGER IDIV,IINT
907C REAL
908 my_real
909 . r0,r1,s0,s,dtol,prec
910C-----------------------------------------------
911C L o c a l V a r i a b l e s
912C-----------------------------------------------
913 INTEGER I
914 my_real
915 . s1,dr,r
916C-----------------------------------------------
917C R1=ABS(E1/E0)
918C IF (R1<=ONE) THEN
919 dr=abs(r0-r1)
920 r=min(r1/prec,dr)
921c R=MIN(R1,DR*PREC)
922 IF ((r<dtol.OR.idiv<-nls_lim)
923 . .OR.(r1>one.AND.iint>0.AND.irefi<2)) THEN
924 IF (imconv==-1)imconv=0
925c IDIV=0
926 idiv=idiv-1
927 RETURN
928 ELSEIF (r1>r0.AND.imconv==-1) THEN
929 idiv=-nls_lim-2
930 s=s0
931 RETURN
932 ENDIF
933 imconv=-1
934 s1=s
935 s=abs(s0*r1-s1*r0)/dr
936 IF (r1>one) THEN
937C----------diverge---------------------------
938 IF (s>one) s=half*s1/r1
939 s=max(em01,s)
940C IDIV=IDIV+1
941 ELSE
942C------extrapolation-------
943 IF (s>max(s1,s0).OR.s<min(s1,s0)) THEN
944 s=min(three_half*s1,s)
945 s=min(five,s)
946 ENDIF
947 ENDIF
948 s0=s1
949 r0=r1
950 idiv=idiv-1
951C
952 RETURN

◆ line_s1()

subroutine line_s1 ( e1,
s,
ep,
sp,
en,
sn,
integer idiv,
dtol,
integer icont,
integer icont0,
integer iriks )

Definition at line 1001 of file nl_solv.F.

1002C-----------------------------------------------
1003C I m p l i c i t T y p e s
1004C-----------------------------------------------
1005#include "implicit_f.inc"
1006C-----------------------------------------------
1007C C o m m o n B l o c k s
1008C-----------------------------------------------
1009#include "impl1_c.inc"
1010C-----------------------------------------------
1011C D u m m y A r g u m e n t s
1012C-----------------------------------------------
1013 INTEGER IDIV,ICONT,ICONT0,IRIKS
1014C REAL
1015 my_real
1016 . e1,s,ep,sp,en,sn,dtol
1017C-----------------------------------------------
1018C L o c a l V a r i a b l e s
1019C-----------------------------------------------
1020 INTEGER I,IM
1021 my_real
1022 . sold,stmp,s2,tol,eold
1023C-------min e=r.du-----------------E1-> already relative /E0-----------------------
1024C------find sl so that el<0
1025 imconv =-1
1026 eold=zero
1027 im = max(icont,icont0)
1028 IF (e1>zero) THEN
1029 IF (idiv==0) THEN
1030 idiv = idiv + 1
1031 sp = one
1032 ep = e1
1033 s = s + one
1034 IF (iriks>0.OR.im > 0) THEN
1035 s = one
1036 imconv =0
1037 END IF
1038 RETURN
1039 ELSEIF (idiv>0) THEN
1040 sold = s
1041 IF (abs(e1-ep)<dtol.OR.e1>ep.OR.idiv>nls_lim) THEN
1042 IF (e1>ep) s = s - one
1043 imconv =0
1044 ELSE
1045 s = s + one
1046 ENDIF
1047 idiv = idiv + 1
1048 sp = sold
1049 ep = e1
1050 RETURN
1051 ENDIF
1052 sp = s
1053 ep = e1
1054 ELSE
1055 IF (idiv==0) THEN
1056 sp = zero
1057 ep = one
1058 s = half
1059 idiv = idiv - 1
1060 sn = one
1061 en = e1
1062 RETURN
1063 ELSEIF (sp==zero) THEN
1064 IF (abs(e1-en)<dtol.OR.idiv<-nls_lim) THEN
1065 imconv =0
1066 ELSEIF (iriks>0.OR.ismdisp>0) THEN
1067 imconv =0
1068 ELSE
1069c SN = S
1070c EN = E1
1071 s = s*half
1072 ENDIF
1073 idiv = idiv - 1
1074 RETURN
1075 ELSEIF (idiv>0) THEN
1076 idiv = -idiv
1077 ENDIF
1078 eold=en
1079 sn = s
1080 en = e1
1081 ENDIF
1082 idiv = idiv - 1
1083c----------return to
1084 sold = s
1085 im = max(icont,icont0)
1086 IF (im>0) THEN
1087c S2 = (EP-EN)/SN+SN
1088c STMP = HALF*(S2-SQRT(S2*S2-FOUR*EP))
1089 stmp = half*(-sp+sn)
1090 ELSE
1091 stmp = ep*(sn-sp)/(ep-en)
1092 ENDIF
1093 s = sp + stmp
1094 stmp =(sold-s)/sold
1095 stmp = max(abs(stmp),abs(e1-eold))
1096 IF (stmp<dtol.OR.idiv<-nls_lim) THEN
1097 imconv =0
1098 ELSEIF (iriks>0.AND.
1099 . (s>two.OR.s<em01)) THEN
1100 imconv =0
1101 ELSE
1102 imconv =-1
1103 ENDIF
1104C
1105 RETURN

◆ nl_solv()

subroutine nl_solv ( integer nddl,
integer, dimension(*) iddl,
integer, dimension(*) ndof,
integer, dimension(*) ikc,
d,
dr,
integer nnz,
integer, dimension(*) iadk,
integer, dimension(*) jdik,
diag_k,
lt_k,
f,
integer nddli,
integer, dimension(*) iadi,
integer, dimension(*) jdii,
diag_i,
lt_i,
integer, dimension(*) itok,
integer, dimension(*) iadm,
integer, dimension(*) jdim,
diag_m,
lt_m,
r02,
dd,
ddr,
integer itask0,
integer it,
integer itc,
ru0,
rold,
integer idiv,
integer inprint,
integer icprec,
integer istop,
e02,
de0,
eimp,
integer, dimension(*) inloc,
integer nddl0,
ls,
u02,
gap,
integer, dimension(*) itab,
integer, dimension(*) fr_elem,
integer, dimension(*) iad_elem,
integer, dimension(*) w_ddl,
a,
ar,
v,
ms,
x,
integer, dimension(*) ipari,
type(intbuf_struct_), dimension(*) intbuf_tab,
integer, dimension(*) num_imp,
integer, dimension(*) ns_imp,
integer, dimension(*) ne_imp,
integer nsrem,
integer nsl,
integer icont,
type(prgraph), dimension(*) graphe,
fac_k,
integer, dimension(*) ipiv_k,
integer nk,
integer nmonv,
integer, dimension(*) imonv,
integer, dimension(*) monvol,
type (surf_), dimension(nsurf) igrsurf,
integer, dimension(*) fr_mv,
volmon,
integer, dimension(*) ibfv,
skew,
xframe,
integer mumps_par,
integer, dimension(*) cddlp,
integer, dimension(*) ind_imp,
integer nbintc,
integer, dimension(*) intlist,
integer, dimension(*) newfront,
integer, dimension(*) isendto,
integer, dimension(*) irecvfrom,
integer, dimension(*) irbe3,
integer, dimension(*) lrbe3,
integer ndiv,
integer icont0,
integer isign,
fext,
dg,
dgr,
dg0,
dgr0,
rfext,
ls1,
integer nodft,
integer nodlt,
integer, dimension(*) irbe2,
integer, dimension(*) lrbe2,
integer idiv0,
relres,
character*1 anew_stif )

Definition at line 53 of file nl_solv.F.

74C-----------------------------------------------
75C M o d u l e s
76C-----------------------------------------------
77 USE dsgraph_mod
78 USE intbufdef_mod
79 USE groupdef_mod
80C-----------------------------------------------
81C I m p l i c i t T y p e s
82C-----------------------------------------------
83#include "implicit_f.inc"
84C-----------------------------------------------
85C C o m m o n B l o c k s
86C-----------------------------------------------
87#if defined(MUMPS5)
88#include "dmumps_struc.h"
89#endif
90#include "com01_c.inc"
91#include "com04_c.inc"
92#include "com08_c.inc"
93#include "impl1_c.inc"
94#include "impl2_c.inc"
95#include "units_c.inc"
96#include "task_c.inc"
97C-----------------------------------------------
98C D u m m y A r g u m e n t s
99C-----------------------------------------------
100C----------resol [K]{X}={F}------
101 INTEGER NDDL ,NNZ ,IADK(*),JDIK(*),IADM(*),JDIM(*),ITASK0,
102 . NDDLI ,IADI(*),JDII(*),ITOK(*),ICPREC,INLOC(*),
103 . NDOF(*),IDDL(*),IKC(*),IDIV ,INPRINT,ISTOP,
104 . IT,ITC,NDDL0,ITAB(*),FR_ELEM(*),IAD_ELEM(*),W_DDL(*)
105 INTEGER NE_IMP(*),NSREM ,NSL,ICONT,ISIGN,IMP1,
106 . IPARI(*) ,NUM_IMP(*),NS_IMP(*), IPIV_K(*), NK
107 INTEGER NMONV,IMONV(*),MONVOL(*),FR_MV(*),
108 . IBFV(*),IND_IMP(*),IRBE3(*) ,LRBE3(*),NDIV ,ICONT0,
109 . IRBE2(*),LRBE2(*)
110 INTEGER NEWFRONT(*),NBINTC,INTLIST(*),ISENDTO(*),IRECVFROM(*),
111 . NODFT ,NODLT,IDIV0
112C REAL
113 my_real
114 . diag_k(*),lt_k(*),f(*),r02,diag_m(*),lt_m(*),
115 . diag_i(*),lt_i(*),d(3,*),dr(3,*),dd(3,*),ddr(3,*),
116 . ru0,rold,e02 ,eimp,lstol,de0,ls(*),u02,gap,rbid
117 my_real
118 . a(3,*),ar(3,*),v(3,*),x(3,*),ms(*),fac_k(*),
119 . volmon(*),skew(*),xframe(*),fext(*),dg(3,*),dgr(3,*),
120 . dg0(3,*),dgr0(3,*),ls1(*) ,rfext, relres(*)
121C-----LS(1):R0,LS(2):S0,LS(3):S1 ; LS1(1): EN, LS1(2):SN;LS1(3):DE
122 TYPE(PRGRAPH) :: GRAPHE(*)
123 INTEGER CDDLP(*)
124 TYPE(INTBUF_STRUCT_) INTBUF_TAB(*)
125 TYPE (SURF_) , DIMENSION(NSURF) :: IGRSURF
126C
127#ifdef MUMPS5
128 TYPE(DMUMPS_STRUC) MUMPS_PAR
129#else
130 ! Fake declaration as DMUMPS_STRUC is shipped with MUMPS
131 INTEGER MUMPS_PAR
132#endif
133C
134 character*1 anew_stif
135#ifdef MUMPS5
136C-----------------------------------------------
137C L o c a l V a r i a b l e s
138C-----------------------------------------------
139C----------resol [K]{X}={F}------
140 INTEGER I ,J ,ND,IP,ITMP,ICOV,ICALU,IER,ISETK0,IRIKS,
141 . F_DDL ,L_DDL
142C REAL
143 my_real
144 . r2,u2,du2,tol,rr,ru,temp,fr,re,de,r1,nl_tol,r01,rr1,
145 . gapn,lamda,fa,uc2,um1,um2,pmax
146C-----------------------------------------------
147C------------------------------------------
148c insolv=1 => (elast) Newton modified
149c insolv=2,3 => BFGS
150c insolv=4,5 => elastoplastic
151C------------------------------------------
152C------------------------------------------
153C------------------------------------------
154C [K]:matrice de rigidite [M]:matrice de preconditionnement
155c F_DDL=1+ITASK*NDDL/NTHREAD
156c L_DDL=(ITASK+1)*NDDL/NTHREAD
157C
158 tol =l_tol
159 nl_tol=n_tol
160C
161 lstol=ls_tol
162 ip=iabs(inprint)
163 de=0
164 nd=3*numnod
165 ru=ru0
166 isetk0=isetk+icont
167 IF (nitol==3.OR.nitol==13.OR.nitol==23
168 . .OR.nitol==123) THEN
169 icalu=1
170 ELSE
171 icalu=0
172 ENDIF
173 iriks=0
174 IF (ncycle>1.AND.idtc==3) THEN
175 IF (ilast==0) iriks=1
176 ENDIF
177C----------------------
178c CALL MY_BARRIER
179C---------------------
180 isetk=0
181 IF (icont>0) gapn=gap*zep9
182 IF (imconv>=0) icont0=icont
183C--------idtc=3 :modified Riks arc_length-----
184 IF (it==0) THEN
185 r01=sqrt(r02)
186C--------ROLD=R2(IT=0)---------------------------------------
187Ctmp IF (ISPRB==1) R01=MIN(R01,SQRT(ROLD)/R01)
188 re=one
189C--------cas particulier---------------------------------------
190 IF (r01*nl_tol<=em12.AND.nl_tol>em10.AND.idiv/=10) THEN
191 IF (irefi/=5.OR.icont==0) THEN
192 IF(inprint/=0)THEN
193 WRITE(iout,*)
194c WRITE(IOUT,1008)IT,RE,R01,RE
195 WRITE(iout,1108)it,anew_stif,re,r01,re
196 WRITE(iout,*)
197 IF(inprint<0)THEN
198 WRITE(istdo,*)
199c WRITE(ISTDO,1008)IT,RE,R01,RE
200 WRITE(istdo,1108)it,anew_stif,re,r01,re
201 WRITE(istdo,*)
202 ENDIF
203 ENDIF
204 relres(1) = re
205 relres(2) = r01
206 relres(3) = re
207 idiv = -2
208 isetk=1
209 RETURN
210 ELSE
211 f(1)=max(em20,f(1))
212 r01=max(r01,em20)
213 ENDIF
214 ENDIF
215 IF (insolv==4) icov = imconv
216 imconv=0
217 rr=rold
218 ndiv=0
219C---------
220 IF (insolv==2.OR.insolv==3) CALL bfgs_0
221C---------
222 IF (ncycle==1) isign = 1
223 IF (iriks>0) THEN
224 CALL lin_solv(nddl ,iddl ,ndof ,ikc ,dg ,
225 1 dgr ,tol ,nnz ,iadk ,jdik ,
226 2 diag_k,lt_k ,nddli ,iadi ,jdii ,
227 3 diag_i,lt_i ,itok ,iadm ,jdim ,
228 4 diag_m,lt_m ,fext ,de ,inloc ,
229 5 fr_elem,iad_elem,w_ddl,itask0,icprec,
230 6 istop ,a ,ar ,v ,
231 7 ms ,x ,ipari ,intbuf_tab ,
232 8 num_imp,ns_imp,ne_imp,nsrem ,nsl ,
233 9 itc ,graphe, itab, fac_k, ipiv_k,
234 a nk ,nmonv,imonv ,monvol,igrsurf,
235 b fr_mv ,volmon,ibfv ,skew ,
236 c xframe ,mumps_par,cddlp,ind_imp,rbid,
237 e irbe3 ,lrbe3 ,irbe2 ,lrbe2 )
238
239C---------
240 CALL produt_uhp(nddl0 ,nddl ,iddl ,ndof ,ikc ,
241 . dg ,dgr ,u2 ,w_ddl )
242C-------distinquer converge or diverge---
243citask0 IF (ITASK == 0) THEN
244 IF (idiv/=-2) THEN
245 CALL cp_real(nd,dd,dg0)
246 IF (iroddl/=0) CALL cp_real(nd,ddr,dgr0)
247 ELSE
248 dla_riks = isign*dt2
249 ENDIF
250citask0 END IF !(ITASK == 0) THEN
251 CALL produt_uhp2(nddl0 ,nddl ,iddl ,ndof ,ikc ,
252 . dg ,dgr ,dg0 ,dgr0 ,uc2 ,
253 . w_ddl )
254C----------------------
255c CALL MY_BARRIER
256C---------------------
257citask0 IF (ITASK == 0) THEN
258 CALL get_max(numnod,dg0,temp,i)
259 um1=temp
260 CALL get_max(numnod,dgr0,temp,i)
261 um1=max(um1,temp,em20)
262 CALL get_max(numnod,dg,temp,i)
263 um2=temp
264 CALL get_max(numnod,dgr,temp,i)
265 um2=max(um2,temp,em20)
266 uc2=uc2/um1/um2
267C LAMDA = ALEN/SQRT(U2)
268 IF ((uc2+dla_riks/rfext)<zero) THEN
269 isign = -1
270 ELSE
271 isign = 1
272 ENDIF
273c print *,'UC2,DLA_RIKS,ISIGN=',UC2,DLA_RIKS,ISIGN
274 IF (impdeb>0) THEN
275 IF (ncycle>=ndeb0.AND.ncycle<=ndeb1)
276 . WRITE(iout,1025)uc2,dla_riks,isign
277 ENDIF
278 IF (isign<0) THEN
279 dla_riks = 2*isign*dt2
280 ELSE
281 dla_riks = zero
282 ENDIF
283C---------
284citask0 END IF !(ITASK == 0) THEN
285
286C----------------------
287c CALL MY_BARRIER
288C---------------------
289 lamda = isign*dt2
290 IF (idiv/=-2.AND.icont==0)THEN
291 CALL produt_vmhp(nddl0 ,nddl ,iddl ,ndof ,ikc ,
292 . dd ,ddr ,f ,de ,w_ddl )
293 CALL integrator2_hp(iddl,ndof,ikc,d,dr,dd,ddr)
294 ELSE
295 CALL frac_dd_hp(iddl,ndof,ikc,d,dr,dg,dgr,lamda)
296 CALL produt_vmhp(nddl0 ,nddl ,iddl ,ndof ,ikc ,
297 . d ,dr ,f ,de ,w_ddl )
298 ENDIF
299 ELSE
300
301 IF ((n_lim /= 1.OR.isprb /= 0).AND.
302 . (ncycle>1.AND.idiv/=-2.AND.icont==0.AND.ikt==0
303 . .OR.(ncycle==1.AND.idiv==10)))THEN
304 IF (idtc==3.AND.ilast==1) THEN
305 lamda = dla_riks
306 dla_riks = zero
307 CALL frac_d_hp(iddl,ndof,ikc,dd,ddr,lamda)
308 ENDIF
309 CALL produt_vmhp(nddl0 ,nddl ,iddl ,ndof ,ikc ,
310 . dd ,ddr ,f ,de ,w_ddl )
311 ELSE
312 IF (insolv==4.AND.icov>=0) imconv=-1
313
314 CALL lin_solv(nddl ,iddl ,ndof ,ikc ,dd ,
315 1 ddr ,tol ,nnz ,iadk ,jdik ,
316 2 diag_k,lt_k ,nddli ,iadi ,jdii ,
317 3 diag_i,lt_i ,itok ,iadm ,jdim ,
318 4 diag_m,lt_m ,f ,de ,inloc ,
319 5 fr_elem,iad_elem,w_ddl,itask0,icprec,
320 6 istop ,a ,ar ,v ,
321 7 ms ,x ,ipari ,intbuf_tab ,
322 8 num_imp,ns_imp,ne_imp,nsrem ,nsl ,
323 9 itc ,graphe, itab, fac_k, ipiv_k,
324 a nk ,nmonv,imonv ,monvol,igrsurf,
325 b fr_mv ,volmon,ibfv ,skew ,
326 c xframe ,mumps_par,cddlp,ind_imp,rbid,
327 e irbe3 ,lrbe3 ,irbe2 ,lrbe2 )
328 IF (de<zero.AND.irefi>1.AND.irefi<5) imconv=-2
329 IF (icont>0.AND.idsgap>0) THEN
330 CALL get_max(numnod,dd,pmax,i)
331C----------------------
332c CALL MY_BARRIER
333C---------------------
334 IF (pmax>gapn) THEN
335 temp=gapn/pmax
336 CALL frac_d_hp(iddl,ndof,ikc,dd,ddr,temp)
337 de = de * temp
338 ENDIF
339 ENDIF
340 ENDIF
341C
342 CALL integrator2_hp(iddl,ndof,ikc,d,dr,dd,ddr)
343 ENDIF !IF (IRIKS>0)
344C----------------------
345c CALL MY_BARRIER
346C---------------------
347citask0 IF (ITASK == 0) THEN
348 IF (imconv>=0) THEN
349 CALL produt_uhp(nddl0 ,nddl ,iddl ,ndof ,ikc ,
350 . d ,dr ,u2 ,w_ddl )
351 eimp=eimp+de
352 IF (de>ep10.AND.iline==0.AND.idyna==0
353 . .AND.ncycle==1) THEN
354 IF(ispmd==0)THEN
355 WRITE(iout,1030)eimp
356 WRITE(istdo,1030)eimp
357 ENDIF
358 CALL imp_stop(-1)
359 ENDIF
360C--------pour idtc>=2-----
361 u02=u2
362 IF(inprint/=0)THEN
363 ru=sqrt(u2)
364c WRITE(IOUT,*)
365c WRITE(IOUT,1002)IT,RU,R01,EIMP
366 WRITE(iout,1102)
367 WRITE(iout,1003)it,anew_stif,ru,r01,eimp
368 IF(inprint<0)THEN
369c WRITE(ISTDO,*)
370c WRITE(ISTDO,1002)IT,RU,R01,EIMP
371 WRITE(istdo,1102)
372 WRITE(istdo,1003)it,anew_stif,ru,r01,eimp
373 ENDIF
374 ENDIF
375
376 relres(1) = ru
377 relres(2) = r01
378 relres(3) = eimp
379 ELSE
380 WRITE(iout,*)
381 WRITE(iout,1013) de
382 IF(inprint<0)THEN
383 WRITE(istdo,*)
384 WRITE(istdo,1013)de
385 ENDIF
386C
387 ENDIF !IF (IMCONV>=0)
388citask0 END IF !(ITASK == 0) THEN
389C---------IT>0---------
390 ELSE
391C--------R=Fext-Fint correction due to Fext---------
392 IF (iriks>0.AND.dla_riks/=zero) THEN
393 CALL vaxpy_hp(nddl, f ,fext,dla_riks)
394 ENDIF
395 CALL produt_hp(nddl,f,f,w_ddl,r2)
396C----------------------
397c CALL MY_BARRIER
398C---------------------
399 IF (r2>=zero.AND.r2<ep30) THEN
400 rr=sqrt(r2/r02)
401 IF (it==1.AND.isprb==1.AND.idiv/=-2) THEN
402 rr1=sqrt(r2/rold)
403C
404 IF (icont>0.AND.rr<one) rr1=rr
405 ELSE
406 rr1=rr
407 ENDIF
408 re=de0/eimp
409 IF (it==1) THEN
410 ru=one
411 ELSEIF (icalu>0) THEN
412 CALL produt_uhp(nddl0 ,nddl ,iddl ,ndof ,ikc ,
413 . d ,dr ,u2 ,w_ddl )
414 CALL produt_uhp(nddl0 ,nddl ,iddl ,ndof ,ikc ,
415 . dd ,ddr ,du2 ,w_ddl )
416C----------------------
417c CALL MY_BARRIER
418C---------------------
419 ru=ls(3)*sqrt(du2/u2)
420 IF (insolv==5) ru=sqrt(du2/u2)
421 ENDIF
422 CALL crit_ite(it,ru,rr1,re,ndiv,nl_tol)
423 ELSE
424 imconv =-2
425 rr1 =ep30
426 ENDIF
427C----------------------
428c CALL MY_BARRIER
429C---------------------
430 IF (imconv==-3) THEN
431 IF (ncycle==1.OR.idiv==-2) imconv=-2
432 r02=r2
433 ENDIF
434C----------------------
435c CALL MY_BARRIER
436C---------------------
437 IF (imconv==0) THEN
438 idiv=0
439C-------start line-search------------
440 IF (iline_s/=1) THEN
441 ls(1)=one
442 ls(2)=zero
443 ENDIF
444 ls(3)=one
445 IF (ndiver>0) r02=r2
446 ENDIF
447 IF(imconv==1) THEN
448 istop=0
449 idiv=0
450C------For time correction-with Riks----
451 IF (it>1.AND.icalu==0) THEN
452 CALL produt_uhp(nddl0 ,nddl ,iddl ,ndof ,ikc ,
453 . d ,dr ,u2 ,w_ddl )
454 CALL produt_uhp(nddl0 ,nddl ,iddl ,ndof ,ikc ,
455 . dd ,ddr ,du2 ,w_ddl )
456C----------------------
457c CALL MY_BARRIER
458C---------------------
459 ru=ls(3)*sqrt(du2/u2)
460 u02=u2
461 ELSEIF (idtc>=2) THEN
462 CALL produt_uhp(nddl0 ,nddl ,iddl ,ndof ,ikc ,
463 . d ,dr ,u02 ,w_ddl )
464C----------------------
465c CALL MY_BARRIER
466C---------------------
467 ENDIF
468C
469 IF(inprint/=0)THEN
470c WRITE(IOUT,*)
471c WRITE(IOUT,1008)IT,RU,RR,RE
472 WRITE(iout,1108)it,anew_stif,ru,rr,re
473 WRITE(iout,*)
474 IF(inprint<0)THEN
475c WRITE(ISTDO,*)
476c WRITE(ISTDO,1008)IT,RU,RR,RE
477 WRITE(istdo,1108)it,anew_stif,ru,rr,re
478 WRITE(istdo,*)
479 ENDIF
480 ENDIF
481 relres(1) = ru
482 relres(2) = rr
483 relres(3) = re
484 ELSEIF(imconv<=-2) THEN
485C--------diverge----
486C write(iout,*)'diverge',IMCONV,IT ,IDIV
487 IF(inprint/=0)THEN
488 WRITE(iout,*)
489 IF (rr1>one) THEN
490 WRITE(iout,1011)rr1
491 ELSEIF(it>nl_dtn) THEN
492 WRITE(iout,1010)
493 ELSEIF(idtc==3) THEN
494 WRITE(iout,1020)
495 ENDIF
496 WRITE(iout,*)
497 IF(inprint<0)THEN
498 WRITE(istdo,*)
499 IF (rr1>one) THEN
500 WRITE(istdo,1011)rr1
501 ELSEIF(it>nl_dtn) THEN
502 WRITE(istdo,1010)
503 ELSEIF(idtc==3) THEN
504 WRITE(istdo,1020)
505 ENDIF
506 WRITE(istdo,*)
507 ENDIF
508 IF(imconv==-2) THEN
509 WRITE(iout,1012)
510 IF(inprint<0)WRITE(istdo,1012)
511 ENDIF
512 ENDIF
513 relres(1) = ru
514 relres(2) = rr1
515 relres(3) = re
516 ELSE
517C--------CAS IMCONV=0,-1---add line-search for constant dt-
518c IF (IDTC>0) THEN
519 IF (iline_s>0.AND.isign>=0.AND.irwall==0) THEN
520 fr=ls(3)
521 IF (nitol/=2.AND.nitol/=4.OR.iline_s==1) THEN
522C----------DD,DDR are modified in some cases
523 IF (ls1(3)/=zero ) THEN
524 de = ls1(3)
525 ELSE
526 CALL produt_vmhp(nddl0 ,nddl ,iddl ,ndof ,ikc ,
527 . dd ,ddr ,f ,de ,w_ddl )
528 END IF
529C----------------------
530c CALL MY_BARRIER
531C---------------------
532 r1=de/de0
533 IF (iline_s==3) THEN
534 r1=abs(de/de0)
535 IF (rr>one) r1=one
536 IF (rr>one.AND.irefi>=2) r1=ls(1)
537 temp=ep02
538 ENDIF
539 ELSE
540 r1=rr/rold
541 temp=ep03
542 ENDIF
543 icov = imconv
544C------------
545citask0 IF (ITASK == 0) THEN
546C------------
547 IF (iline_s==3) THEN
548 CALL line_s(ls(1),ls(2),r1,fr,lstol,idiv,icont,temp)
549 ELSEIF (iline_s==1) THEN
550 CALL line_s1(r1,fr,ls(1),ls(2),ls1(1),ls1(2),idiv,lstol,
551 . icont,icont0,iriks)
552 ELSEIF (iline_s==2) THEN
553 CALL line_s(ls(1),ls(2),r1,fr,lstol,idiv,icont,temp)
554 ENDIF
555 IF (impdeb>0.AND.ispmd==0) THEN
556 IF (ncycle>=ndeb0.AND.ncycle<=ndeb1) then
557 WRITE(iout,*)'R1,FR,IDIV=',r1,fr,idiv
558C WRITE(IOUT,*)'DE,DE0=',DE,DE0
559 end if
560 ENDIF
561C------------
562citask0 END IF !(ITASK == 0) THEN
563C----------------------
564c CALL MY_BARRIER
565C---------------------
566C
567 IF ((insolv==2.OR.insolv==3).AND.fr>one.AND.ikt==0)
568 . fr = min(fr,two)
569 ENDIF
570C IF (IABS(INSOLV)==5) THEN
571C----------suprimed----------
572C ELSE
573 IF (imconv==0) THEN
574C----------schema normale---------
575 IF (it>1) THEN
576C------------
577 CALL produt_uhp(nddl0 ,nddl ,iddl ,ndof ,ikc ,
578 . d ,dr ,u2 ,w_ddl )
579 CALL produt_uhp(nddl0 ,nddl ,iddl ,ndof ,ikc ,
580 . dd ,ddr ,du2 ,w_ddl )
581C----------------------
582c CALL MY_BARRIER
583C---------------------
584 ru=ls(3)*sqrt(du2/u2)
585 ENDIF
586C------------
587citask0 IF (ITASK == 0) THEN
588C------------
589 IF (insolv==2.OR.insolv==3) CALL bfgs_ls(ls(3))
590 IF (rr>=rold) THEN
591 ndiv = ndiv +1
592 ELSE
593 ndiv = 0
594 ENDIF
595 IF (inprint/=0)THEN
596 IF(mod(it,ip)==0)THEN
597c WRITE(IOUT,*)
598 WRITE(iout,1003)it,anew_stif,ru,rr,re
599 IF(inprint<0)THEN
600c WRITE(ISTDO,*)
601 WRITE(istdo,1003)it,anew_stif,ru,rr,re
602 ENDIF
603 ENDIF
604 ENDIF
605 relres(1) = ru
606 relres(2) = rr
607 relres(3) = re
608 IF(itc>=n_lim.AND.ismdisp==0) THEN
609c IF(INPRINT/=0)THEN
610c WRITE(IOUT,*)
611c WRITE(IOUT,1006)ITC
612c IF(INPRINT<0)THEN
613c WRITE(ISTDO,*)
614c WRITE(ISTDO,1006)ITC
615c ENDIF
616c ENDIF
617 itc=0
618 ENDIF
619C------------
620citask0 END IF !(ITASK == 0) THEN
621C----------------------
622c CALL MY_BARRIER
623C---------------------
624 IF (iriks>0) THEN
625 IF (isetk0==1) THEN
626 CALL lin_solv(nddl ,iddl ,ndof ,ikc ,dg ,
627 1 dgr ,tol ,nnz ,iadk ,jdik ,
628 2 diag_k,lt_k ,nddli ,iadi ,jdii ,
629 3 diag_i,lt_i ,itok ,iadm ,jdim ,
630 4 diag_m,lt_m ,fext ,de ,inloc ,
631 5 fr_elem,iad_elem,w_ddl,itask0,icprec,
632 6 istop ,a ,ar ,v ,
633 7 ms ,x ,ipari ,intbuf_tab ,
634 8 num_imp,ns_imp,ne_imp,nsrem ,nsl ,
635 9 itc ,graphe, itab, fac_k, ipiv_k,
636 a nk ,nmonv,imonv ,monvol,igrsurf,
637 b fr_mv ,volmon,ibfv ,skew ,
638 c xframe ,mumps_par,cddlp,ind_imp,rbid,
639 e irbe3 ,lrbe3 ,irbe2 ,lrbe2 )
640 icprec = 0
641 idsc = 0
642 END IF
643 CALL produt_vmhp(nddl0 ,nddl ,iddl ,ndof ,ikc ,
644 . dg ,dgr ,f ,um1 ,w_ddl )
645 END IF
646C
647 CALL lin_solv(nddl ,iddl ,ndof ,ikc ,dd ,
648 1 ddr ,tol ,nnz ,iadk ,jdik ,
649 2 diag_k,lt_k ,nddli ,iadi ,jdii ,
650 3 diag_i,lt_i ,itok ,iadm ,jdim ,
651 4 diag_m,lt_m ,f ,de ,inloc ,
652 5 fr_elem,iad_elem,w_ddl,itask0,icprec,
653 6 istop ,a ,ar ,v ,
654 7 ms ,x ,ipari ,intbuf_tab ,
655 8 num_imp,ns_imp,ne_imp,nsrem ,nsl ,
656 9 itc ,graphe, itab, fac_k, ipiv_k,
657 a nk ,nmonv,imonv ,monvol,igrsurf,
658 b fr_mv ,volmon,ibfv ,skew ,
659 c xframe ,mumps_par,cddlp,ind_imp,rbid,
660 d irbe3 ,lrbe3 ,irbe2 ,lrbe2 )
661 IF (icont>0.AND.idsgap>0) THEN
662 CALL get_max(numnod,dd,pmax,i)
663C----------------------
664c CALL MY_BARRIER
665C---------------------
666 IF (pmax>gapn) THEN
667 temp=gapn/pmax
668 CALL frac_d_hp(iddl,ndof,ikc,dd,ddr,temp)
669 de = de * temp
670 ENDIF
671 ENDIF
672C----------------------
673c CALL MY_BARRIER
674C---------------------
675C---------Riks type arc length------
676 IF (iriks>0) THEN
677 lamda = dla_riks
678 IF (ial_m==1) THEN
679 CALL al_constraint1_hp(nddl0 ,nddl ,iddl ,ndof ,ikc ,
680 . dd ,ddr ,dg ,dgr ,d ,
681 . dr ,w_ddl ,alen ,lamda ,scal_riks,
682 . ier )
683C----------------------
684c CALL MY_BARRIER
685C---------------------
686 IF (ier==1) THEN
687 imconv=-2
688 IF(inprint/=0)THEN
689 WRITE(iout,*)
690 WRITE(iout,1020)
691 WRITE(iout,*)
692 IF(inprint<0)THEN
693 WRITE(istdo,*)
694 WRITE(istdo,1020)
695 WRITE(istdo,*)
696 ENDIF
697 ENDIF
698 END IF !(IER==1) THEN
699 ELSEIF (ial_m==2) THEN
700 CALL al_constraint2_hp(nddl0 ,nddl ,iddl ,ndof ,ikc ,
701 . dd ,ddr ,dg ,dgr ,d ,
702 . dr ,w_ddl ,alen ,lamda ,scal_riks)
703C----------------------
704c CALL MY_BARRIER
705C---------------------
706 END IF
707 IF (imconv>=0) THEN
708 CALL frac_dd_hp(iddl,ndof,ikc,d ,dr ,dg,dgr,lamda)
709citask0 IF (ITASK == 0) THEN
710 dla_riks = dla_riks + lamda
711 eimp=eimp+dla_riks*um1
712citask0 END IF !(ITASK == 0) THEN
713C----------------------
714c CALL MY_BARRIER
715C---------------------
716 ENDIF
717 END IF !IF (IRIKS>0) THEN
718C
719 CALL integrator2_hp(iddl,ndof,ikc,d,dr,dd,ddr)
720 eimp=eimp+de
721c if (de<0) print *,'****negative de=',de
722C
723 ELSEIF(imconv==-1) THEN
724C----------line-search---------
725 temp=-ls(3)+fr
726 CALL frac_dd_hp(iddl,ndof,ikc,d,dr,dd,ddr,temp)
727C----------------------
728c CALL MY_BARRIER
729C---------------------
730 ls(3) = fr
731 END IF !IF (IMCONV==0) THEN
732C
733 END IF !IF(IMCONV==1) THEN
734C-----fin if IT=0
735 ENDIF
736 IF (imconv>=0) THEN
737 ru0=ru
738 de0=de
739 IF (abs(de0)<em20) de0=em20
740 rold=max(em20,rr)
741 ENDIF
742C-----REFORMING [K]----
743 IF (ismdisp>0) THEN
744 IF (itc>=n_lim.AND.imconv>=0) THEN
745 itc = 0
746 isetk=1
747 END IF
748 IF (itc==0.AND.imconv>=0) THEN
749 IF (ncycle==1.OR.n_lim==1) isetk=1
750 IF (isolv==5.OR.isolv==6) isetk=1
751 END IF
752 IF (irwall >0 .AND.imconv == 1) isetk=1
753 IF ((itc==0.OR.imconv==1).AND.ikt>0) isetk=1
754 ELSE
755C----condition IMCONV==1 reforming [K] for nothing (no resolution)-
756 IF (itc==0.OR.imconv==1) isetk=1
757 END IF
758 IF ((idtc==3.OR.ikt>0).AND.imconv<=-2) isetk=1
759C-----if still diverge, try w/o [K] reforming
760 IF (imconv<=-2.AND.rr1>ep20.AND.idiv0>=0) isetk=1
761 IF ((istop==1.OR.istop==2).AND.(isolv==5.OR.isolv==6))THEN
762 istop = 0
763 isetk = 1
764 imconv=-2
765 ENDIF
766 IF (imconv<=-2.AND.isprb>0.AND.rr1>ep04) isetk=1
767! sb
768 IF(imconv<= -2.AND.n_lim == 1 .AND. isprb == 0) isetk = 1
769
770C----------------------
771c CALL MY_BARRIER
772C---------------------
773Ctmp IF (ITC==0.OR.IMCONV==1.OR.IMCONV<=-2) ISETK=1
774c 1002 FORMAT(5X,'NL_ITERATION=',I5,1X,' INITIAL RESIDUAL NORM=',3E11.4)
775 1102 FORMAT(3x,78('-')/
776 . 12x,'Stif. Mat.'/
777 . 6x,'Iter',2x, ' reformed ',5x,'|du|/|u|',3x,'|r|/|r0|',3x,'|dE|/|E|',1x,'Conv.stat.'/
778 . 3x,78('-'))
779c 1003 FORMAT(5X,'NL_ITERATION=',I5,1X,' RELATIVE RESIDUAL NORM=',3E11.4)
780 1003 FORMAT(5x,i5,6x,a1,5x,3(1x,1pe10.3))
781 1005 FORMAT(3x,'TOLERANCE FOR LINEAR ITERATIVE SOLVER :',2x,e11.4)
782 1006 FORMAT(3x,'--STIFFNESS MATRIX IS RESET AFTER ',i4,
783 . ' ITERATIONS--')
784 1007 FORMAT(3x,'--ITERATION DIVERGE with R=',e11.4,2x,
785 . 'STIFFNESS MATRIX WILL BE RESET')
786c 1008 FORMAT(3X,'* CONVERGED WITH ',I4,1X,'ITERATIONS,',
787c . ' |du|/|u|,|r|/|r0|,|dE|/|E|=',3E11.4)
788 1108 FORMAT(5x,i5,6x,a1,5x,3(1x,1pe10.3),5x, 'C')
789 1009 FORMAT(3x,'--ITERATION DIVERGE with RELATIVE R=',e11.4/,
790 . 3x,'RESET ITERATION WITH STABILIZATION BY SCALING= ',e11.4)
791 1010 FORMAT(3x,'--ITERATION DIVERGE with MAX_ITER REACHED--')
792 1011 FORMAT(3x,'ITERATION DIVERGE with RELATIVE R=',e11.4)
793 1012 FORMAT(3x,'--RESET ITERATION WITH NEW TIMESTEP--')
794 1013 FORMAT(3x,'ITERATION DIVERGE with NEGATIVE ENERGY DE=',e11.4)
795 1020 FORMAT(3x,'--ITERATION DIVERGE with RIKS Method--')
796 1025 FORMAT(3x,'UC2,DLA_RIKS,ISIGN=',2e11.4,i4)
797 1030 FORMAT(3x,'CHECK CONSTRAINT CONDITIONS, TOO LARGE ENERGY VALUE=',
798 . 2x,e11.4)
799 RETURN
800#endif
bfgs_0
subroutine bfgs_0
Definition imp_bfgs.F:81
bfgs_ls
subroutine bfgs_ls(ls)
Definition imp_bfgs.F:109
imp_stop
subroutine imp_stop(istop)
Definition imp_solv.F:1992
frac_d_hp
subroutine frac_d_hp(iddl, ndof, ikc, d, dr, fac)
Definition integrator.F:565
lin_solv
subroutine lin_solv(nddl, iddl, ndof, ikc, d, dr, tol, nnz, iadk, jdik, diag_k, lt_k, nddli, iadi, jdii, diag_i, lt_i, itok, iadm, jdim, diag_m, lt_m, f, f_u, inloc, fr_elem, iad_elem, w_ddl, itask, icprec, istop, a, ar, ve, ms, xe, ipari, intbuf_tab, num_imp, ns_imp, ne_imp, nsrem, nsl, it, graphe, itab, fac_k, ipiv_k, nk, nmonv, imonv, monvol, igrsurf, fr_mv, volmon, ibfv, skew, xframe, mumps_par, cddlp, ind_imp, xi_c, irbe3, lrbe3, irbe2, lrbe2)
Definition lin_solv.F:74
dsgraph_mod
Definition dsgraph_mod.F:31
groupdef_mod
Definition groupdef_mod.F:662
al_constraint2_hp
subroutine al_constraint2_hp(nddl0, nddl, iddl, ndof, ikc, dd, ddr, dg, dgr, di, dir, w_ddl, l_a, lamda, sw2)
Definition nl_solv.F:1541
crit_ite
subroutine crit_ite(it, ur, rr, er, ndiv, tol)
Definition nl_solv.F:809
al_constraint1_hp
subroutine al_constraint1_hp(nddl0, nddl, iddl, ndof, ikc, dd, ddr, dg, dgr, di, dir, w_ddl, l_a, lamda, sw2, ier)
Definition nl_solv.F:1463
get_max
subroutine get_max(n, d, vmax, i)
Definition nl_solv.F:962
line_s1
subroutine line_s1(e1, s, ep, sp, en, sn, idiv, dtol, icont, icont0, iriks)
Definition nl_solv.F:1002
line_s
subroutine line_s(r0, s0, r1, s, dtol, idiv, iint, prec)
Definition nl_solv.F:895
cp_real
subroutine cp_real(n, x, xc)
Definition produt_v.F:871
produt_hp
subroutine produt_hp(nddl, x, y, w, r)
Definition produt_v.F:3252
vaxpy_hp
subroutine vaxpy_hp(n, v, y, s)
Definition produt_v.F:3580
produt_vmhp
subroutine produt_vmhp(nddl0, nddl, iddl, ndof, ikc, dd, ddr, y, r, w_imp)
Definition produt_v.F:3321