Télécharger aft.procedur

Retour à la liste

Numérotation des lignes :

  1. * AFT PROCEDUR BP208322 16/09/01 21:15:01 9010
  2. *
  3. ************************************************************************
  4. *
  5. * OBJET : Calculer la force non-lineaire FNL et la derivee
  6. * KNL = -dFNL/dX en basculant entre domaine temporel et
  7. * frequentiel
  8. * TFR : transformee de fourier rapide d'un listreel/listchpo
  9. *
  10. * ENTREE :
  11. * TABAFT = TABLE
  12. * . 'WTABLE' . 'DEPLACEMENTS' = deplacement frequentiel actuel
  13. * . 'N_HARMONIQUE' = nombre d'harmonique
  14. * . 'COEFF_FNL' = coefficient de la force nl
  15. * . 'HARM_FORCE' = comp force par harmonique
  16. * . 'HARM_COMP_DEP' = comp dep par composante
  17. * . 'HARM_COMP_FOR' = comp force par composante
  18. * . 'POINT_FNL' = point ou applique la force nl
  19. * . 'COMP_FNL' = composante de la force nl
  20. * LMOT1 = LISTMOTS de mots-clés définissant les actions a effectuer
  21. * a choisir parmi { RESI, DRES }
  22. *
  23. * PROCEDURE APPELE: CHARMECA
  24. *
  25. * SORTIE :
  26. * TNL1 = TABLE
  27. * . 'ADDI_SECOND' = force non-lineaire frequentiel
  28. * . 'ADDI_MATRICE' = derivee de la force nl frequentiel
  29. *
  30. * CREATION : L. Xie, 2014
  31. * MODIFS : BP, 2016
  32. *
  33. * TODO (BP, 2016) : a terme, cette procedure devra etre transformee en
  34. * un operateur
  35. *
  36. ************************************************************************
  37. DEBPROC AFT TABAFT*'TABLE' time*'FLOTTANT' LMOTAFT/'LISTMOTS' ;
  38.  
  39.  
  40. ************************************************************************
  41. * VERIFICATION DES DONNEES D'ENTREE + PREPARATIFS *
  42. ************************************************************************
  43.  
  44. *coefficient de la force nl
  45. *
  46. * mess (chai 'AVANT AFT temps passé : ' (temp 'NOEC'));
  47. IFCHPO = TABAFT . 'CALC_CHPO';
  48. SI (EXIS TABAFT . 'POINT_FNL');
  49. Pnl = TABAFT . 'POINT_FNL';
  50. SINON;
  51. MESS 'AFT: IL MANQUE LA GEOMETRIE NL DANS POINT_FNL'; ERRE 21;
  52. FINSI;
  53. ttype = TYPE Pnl;
  54. SI (EGA ttype 'POINT');
  55. nPnl = 1;
  56. FINSI;
  57. SI (EGA ttype 'MAILLAGE');
  58. ltype = ELEM Pnl 'TYPE';
  59. si ( ((dime ltype) neg 1) OU (neg (extr ltype 1) 'POI1') );
  60. mess 'AFT: POINT_FNL doit etre un POINT ou un MAILLAGE de POI1';
  61. ERRE 21;
  62. finsi;
  63. nPnl = NBEL Pnl;
  64. FINSI;
  65.  
  66. SI (EXIS TABAFT . 'COMP_FNL');
  67. COMP_FNL = TABAFT . 'COMP_FNL';
  68. n_CNL = DIME COMP_FNL;
  69. SINON;
  70. MESS 'AFT: IL MANQUE LA COMPOSANTE NL dans COMP_FNL'; ERRE 21;
  71. FINSI;
  72.  
  73. SI (EXIS TABAFT 'WTABLE');
  74. WTAB = TABAFT . 'WTABLE';
  75. SINON;
  76. MESS 'IL MANQUE LA WTABLE'; ERRE 21;
  77. FINSI;
  78.  
  79. *DEPLACEMENTS FREQUENTIELS
  80. DEPTOT = WTAB . 'DEPLACEMENTS';
  81. * mess 'aft: t='time ' et DEPTOT='; list DEPTOT;
  82. *
  83. SI (EXIS TABAFT 'N_HARMONIQUE');
  84. nhbm = TABAFT . 'N_HARMONIQUE';
  85. SINON;
  86. MESS 'IL MANQUE LE NOMBRE D HARMONIQUE'; ERRE 21;
  87. FINSI;
  88. *
  89. SI (EXIS TABAFT . 'COMPOSANTES' 'DEPLACEMENT');
  90. NOMCOM = TABAFT . 'COMPOSANTES' . 'DEPLACEMENT';
  91. ncom = DIME NOMCOM ;
  92. SINON;
  93. MESS 'IL MANQUE LES NOMS DE COMPOSANTES'; ERRE 21;
  94. FINSI;
  95.  
  96. SI (EXIS LMOTAFT);
  97. * IFFNL = EXIS LMOTAFT 'RESI';
  98. * IFKNL = EXIS LMOTAFT 'DRES';
  99. IFFNL = (EXIS LMOTAFT 'RESI') ou (EXIS LMOTAFT 'DRES');
  100. IFKNL = EXIS LMOTAFT 'RIGI';
  101. SINON;
  102. IFFNL = VRAI;
  103. IFKNL = FAUX;
  104. FINSI;
  105.  
  106. * faut-il calculer les vitesses ?
  107. FLVITE = TAB1 . 'CALC_VITE';
  108.  
  109. * pour la vitesse, on a besoin de la frequence w
  110. IFREQ = WTAB . 'FREQ_INCONNUE';
  111. si IFREQ; wrads = WTAB . 'FREQ';
  112. sinon; wrads = IPOL TAB1 . 'FREQUENCE' time;
  113. finsi;
  114. * wHz = IPOL TAB1 . 'FREQUENCE' time;
  115. * wrads = 2.*pi * wHz;
  116.  
  117. *liste mots par harmonique
  118. NHBMU = TABAFT . 'COMPOSANTES' . 'HARM_DEPLACEMENT';
  119. NOMU0 = NHBMU . 0;
  120. NHBMF = TABAFT . 'COMPOSANTES'. 'HARM_FORCE';
  121. NOMF0 = NHBMF . 0;
  122. COMPU = TABAFT . 'COMPOSANTES'. 'HARM_COMP_DEP';
  123. COMPF = TABAFT . 'COMPOSANTES'. 'HARM_COMP_FOR';
  124.  
  125.  
  126. ************************************************************************
  127. * PARAMETRES TEMPORELS / FREQUENTIELS *
  128. ************************************************************************
  129.  
  130. * timew = time * 2. * pi;
  131. *bp tdegs = 360. * time ;
  132. tdegs = 360. ;
  133.  
  134. *nombre de point pendant une periode pour faire TFR
  135. OTFR = TABAFT . 'N_PT_TFR';
  136. NPTFR = 2**OTFR;
  137. *bp TDT = 1. / time ;
  138. TDT = 1.;
  139. dtDT = TDT / NPTFR ;
  140.  
  141. * LISTES DES FONCTIONS TRIGONOMETRIQUES :
  142. *
  143. * LISTDT = listreel{ 0. (1./2**OTFR) ... 1. }
  144. * LISDEG = listreel{ 0. (360./2**OTFR) ... 360. }
  145. * pour u(t) :
  146. * LISCOS . k = listreel{ ... cos(k*t_i) ... }
  147. * LISSIN . k = listreel{ ... sin(k*t_i) ... }
  148. * TGAMMA . i = listreel{ 1 ... cos(k*t_i) sin(k*t_i) ... }
  149. * pour v(t)=du(t)/dt :
  150. * LISCOS2 . k = listreel{ ... -k*sin(k*t_i) ... }
  151. * LISSIN2 . k = listreel{ ... k*cos(k*t_i) ... }
  152. * TGAMMA_V . i = listreel{ 0 ... -k*sin(k*t_i) k*cos(k*t_i) ... }
  153.  
  154. * on recupere si les listreels existent deja
  155. SI (EXIS WTAB 'LISTDT');
  156. LISTDT = WTAB . 'LISTDT';
  157. LISDEG = WTAB . 'LISDEG';
  158. LISCOS = WTAB . 'LISCOS';
  159. LISSIN = WTAB . 'LISSIN';
  160. TGAMMA = WTAB . 'TGAMMA';
  161. SI FLVITE;
  162. LISCOS2 = WTAB . 'LISCOS_V';
  163. LISSIN2 = WTAB . 'LISSIN_V';
  164. TGAMMA_V = WTAB . 'TGAMMA_V';
  165. FINSI;
  166.  
  167. * on cree les listreels la premiere fois
  168. SINON;
  169. LISTDT = PROG 0. pas dtDT (TDT) ;
  170. LISDEG = tdegs * LISTDT;
  171. * ---cas LISTREEL et cas calcul de la jacobienne KNL
  172. LISCOS = TABLE;
  173. LISSIN = TABLE;
  174. LISCOS . 0 = COS (0 * LISDEG);
  175. REPE BCOS nhbm; ik = &BCOS;
  176. LISCOS . ik = COS (ik * LISDEG);
  177. LISSIN . ik = SIN (ik * LISDEG);
  178. FIN BCOS;
  179. WTAB . 'LISCOS' = LISCOS;
  180. WTAB . 'LISSIN' = LISSIN;
  181. * et cas calcul de la jacobienne CNL
  182. SI FLVITE;
  183. LISCOS2 = TABLE;
  184. LISSIN2 = TABLE;
  185. LISCOS2 . 0 = SIN (0. * LISDEG);
  186. LISSIN2 . 0 = LISCOS2 . 0;
  187. REPE BCOS nhbm; ik = &BCOS;
  188. LISCOS2 . ik = (-1.* ik) * (LISSIN . ik);
  189. LISSIN2 . ik = ( ik) * (LISCOS . ik);
  190. FIN BCOS;
  191. WTAB . 'LISCOS_V' = LISCOS2;
  192. WTAB . 'LISSIN_V' = LISSIN2;
  193. FINSI;
  194. * --- cas LISTCHPO
  195. TGAMMA = TABLE;
  196. * it = temps; ik = harmonique;
  197. it = 0;
  198. * -BOUCLE SUR LE TEMPS------------------------------
  199. REPE btemps (dime LISTDT); it = it + 1;
  200. xdeg = extr LISDEG it;
  201. * harmonique 0
  202. xt = prog 1.;
  203. ik = 0;
  204. * ---BOUCLE SUR LES HARMONIQUES-----------------------
  205. REPE bouhbm nhbm; ik = ik + 1 ;
  206. xt = xt et (cos (ik*xdeg)) et (sin (ik*xdeg));
  207. FIN bouhbm;
  208. TGAMMA . it = xt;
  209. FIN btemps;
  210. WTAB . 'LISTDT' = LISTDT;
  211. WTAB . 'LISDEG' = LISDEG;
  212. WTAB . 'TGAMMA' = TGAMMA;
  213. * idem pour la VITESSE
  214. SI FLVITE;
  215. TGAMMA_V = TABLE;
  216. it = 0;
  217. REPE btemps (dime LISTDT); it = it + 1;
  218. xdeg = extr LISDEG it;
  219. xt = prog 0.;
  220. ik = 0;
  221. REPE bouhbm nhbm; ik = ik + 1 ;
  222. xt = xt et (-1.*ik*(sin (ik*xdeg))) et (ik*(cos (ik*xdeg)));
  223. FIN bouhbm;
  224. TGAMMA_V . it = xt;
  225. FIN btemps;
  226. WTAB . 'TGAMMA_V' = TGAMMA_V;
  227. FINSI;
  228.  
  229. FINSI;
  230.  
  231.  
  232.  
  233. ************************************************************************
  234. * CALCUL DES DEPLACEMENTS NL TEMPORELLES *
  235. * A PARTIR DES DEPLACEMENTS FREQUENTIELS *
  236. ************************************************************************
  237.  
  238. * calcul via l'appel a AFT1
  239. DEPnl VITnl = AFT1 DEPTOT;
  240.  
  241.  
  242. ************************************************************************
  243. * CALCUL DE LA FORCE NL DANS LE DOMAINE TEMPOREL *
  244. * (APPEL A LA FONCTION UTILISATEUR CHARMECA) *
  245. ************************************************************************
  246.  
  247. * DONNEES D'ENTREE
  248. TABAFT . 'DEP_NL' = DEPnl;
  249. TABAFT . 'VIT_NL' = VITnl;
  250. TABAFT . 'TIME' = TIME;
  251.  
  252. * APPEL A CHARMECA
  253. CHARMECA TABAFT;
  254.  
  255. * DONNEES DE SORTIE
  256.  
  257. * fnl(t) (domaine temporel)
  258. FNLTEMP = TABAFT . 'FNL_T';
  259. typFNL = type FNLTEMP;
  260. IFCHPO2 = ega typFNL 'LISTCHPO';
  261. si (non IFCHPO2);
  262. si (neg typFNL 'TABLE');
  263. mess 'AFT: FNL_T devrait etre de type LISTCHPO ou TABLE!'; erre 21;
  264. finsi;
  265. finsi;
  266.  
  267. * Knl=dfnl/du(t) (domaine temporel)
  268. KNLT = TABAFT . 'KNL_T';
  269. SI FLVITE; CNLT = TABAFT . 'CNL_T'; FINSI;
  270. * c'est toujours une table pour l'instant
  271.  
  272.  
  273.  
  274. ************************************************************************
  275. * CALCUL DE LA FORCE NON-LINEAIRE ET DE SA DERIVEE *
  276. * DANS LE DOMAINE FREQUENTIEL *
  277. ************************************************************************
  278.  
  279. * calcul via l'appel a AFT2
  280. FNLFRE = AFT2 FNLTEMP;
  281.  
  282.  
  283. *=== CALCUL DE KNL et CNL =============================================
  284. *
  285. * TODO (BP) : Il faudrait remplacer les lignes ci-dessous par une
  286. * nouvelle options de TFR ?
  287. *
  288. SI IFKNL;
  289.  
  290. *=== via la TFR de knl et cnl exprimes en temporel ====================
  291. *
  292. *=== CALCUL DE KNL PAR LES LISTREELs ==================================
  293. *
  294. si (neg (dime KNLT) 0);
  295.  
  296. * preparation des objets resultats
  297. * KNLFRE = VIDE 'RIGIDITE';
  298. LISTR1 = prog;
  299.  
  300. *----------------------- boucle sur les lignes de la Jacobienne Knl=-dfnl/du
  301. * (= boucle sur les ddl temporels)
  302. REPE Bligne (DIME KNLT); iligne = &Bligne;
  303.  
  304. KNLT_i = KNLT . iligne;
  305. nKNLT_i = DIME KNLT_i;
  306.  
  307. * harmonique 0 : TFR(Knl=-dfnl/du)
  308. KNLFR = TABLE; KNLFI = TABLE;
  309. * boucle sur les colonnes de la Jacobienne temporelle : dfnl_iligne/du_icnl
  310. REPE Bcnl1 nKNLT_i; icnl = &Bcnl1;
  311. EVKNLT = EVOL MANU LISTDT KNLT_i . icnl;
  312. EVKNLF = TFR OTFR EVKNLT 'REIM' ;
  313. KNLFR . icnl = EXTR EVKNLF 'ORDO' 1;
  314. KNLFI . icnl = EXTR EVKNLF 'ORDO' 2;
  315. FIN Bcnl1;
  316.  
  317. * harmonique ik :
  318. * TFR(Knl * cos(ik wt)) = KNLCR.ik + i KNLCI.ik
  319. * TFR(Knl * sin(ik wt)) = KNLSR.ik + i KNLSI.ik
  320. * tables de listreels
  321. KNLCR = TABLE; KNLCI = TABLE;
  322. KNLSR = TABLE; KNLSI = TABLE;
  323. ik = 0;
  324. REPE Bdtfr nhbm; ik = ik + 1;
  325. DFFCR = TABLE; DFFCI = TABLE;
  326. DFFSR = TABLE; DFFSI = TABLE;
  327. * boucle sur les colonnes de la Jacobienne temporelle : -dfnl_iligne/du_icnl
  328. REPE Bcnl2 nKNLT_i; icnl = &Bcnl2;
  329. KNLTCik = KNLT_i . icnl * LISCOS . ik;
  330. KNLTSik = KNLT_i . icnl * LISSIN . ik;
  331. EVDFTC = EVOL MANU LISTDT KNLTCik;
  332. EVDFTS = EVOL MANU LISTDT KNLTSik;
  333. EVDFFC = TFR OTFR EVDFTC 'REIM';
  334. EVDFFS = TFR OTFR EVDFTS 'REIM';
  335. DFFCR . icnl = EXTR EVDFFC 'ORDO' 1;
  336. DFFCI . icnl = EXTR EVDFFC 'ORDO' 2;
  337. DFFSR . icnl = EXTR EVDFFS 'ORDO' 1;
  338. DFFSI . icnl = EXTR EVDFFS 'ORDO' 2;
  339. FIN Bcnl2;
  340. KNLCR . ik = DFFCR; KNLCI . ik = DFFCI;
  341. KNLSR . ik = DFFSR; KNLSI . ik = DFFSI;
  342. FIN Bdtfr;
  343.  
  344. * on range les valeurs dans LISTR1
  345. * -les n*(2H+1) premieres valeurs
  346. * --> ligne 1 : dF_0/dU = [ dF_0/dU_0 .. dF_0/dU_ik dF_0/dV_ik .. ]
  347. ik = 0;
  348. REPE BML1 nKNLT_i; icnl = &BML1;
  349. LISTR1 = LISTR1 ET ((EXTR KNLFR . icnl 1) );
  350. REPE BMLk nhbm; ik = &BMLk;
  351. LISTR1 = LISTR1 ET ((EXTR KNLCR . ik . icnl 1) );
  352. LISTR1 = LISTR1 ET ((EXTR KNLSR . ik . icnl 1) );
  353. FIN BMLk;
  354. FIN BML1;
  355. *
  356. * -lignes (2,3) , (4,5) ... (cos,sin) ... : dF_ih/dU
  357. REPE BMKk nhbm; ih = &BMKk;
  358. * ligne cos : dF_ih/dU = [ dF_ih/dU_0 .. dF_ih/dU_ik dF_ih/dV_ik .. ]
  359. REPE BMKC1 nKNLT_i; icnl = &BMKC1;
  360. LISTR1 = LISTR1 ET ((EXTR KNLFR . icnl (ih + 1)) * 2 );
  361. REPE BMKCK nhbm; ik = &BMKCK;
  362. LISTR1 = LISTR1 ET
  363. ((EXTR KNLCR . ik . icnl (ih + 1)) * 2 );
  364. LISTR1 = LISTR1 ET
  365. ((EXTR KNLSR . ik . icnl (ih + 1)) * 2 );
  366. FIN BMKCK;
  367. FIN BMKC1;
  368. * ligne sin : dG_ih/dU = [ dG_ih/dU_0 .. dG_ih/dU_ik dG_ih/dV_ik .. ]
  369. REPE BMKS1 nKNLT_i; icnl = &BMKS1;
  370. LISTR1 = LISTR1 ET ((EXTR KNLFI . icnl (ih + 1)) * -2 );
  371. REPE BMKSK nhbm; ik = &BMKSK;
  372. LISTR1 = LISTR1 ET
  373. ((EXTR KNLCI . ik . icnl (ih + 1)) * -2 );
  374. LISTR1 = LISTR1 ET
  375. ((EXTR KNLSI . ik . icnl (ih + 1)) * -2 );
  376. FIN BMKSK;
  377. FIN BMKS1;
  378. FIN BMKk;
  379.  
  380.  
  381. FIN Bligne;
  382. *----------------------------------------- fin de boucle sur les lignes
  383.  
  384. RIDEP = MOTS;
  385. * RIFOR = MOTS;
  386. inl = 0;
  387. REPE Bcfnl1 (nKNLT_i / nPnl); inl = inl + 1;
  388. iCOMNL = EXTR COMP_FNL inl;
  389. RIDEP = RIDEP ET COMPU . iCOMNL ;
  390. * iCOMFNL = CHAI 'F' (EXTR iCOMNL 2);
  391. * RIFOR = RIFOR ET COMPF . iCOMFNL ;
  392. FIN Bcfnl1;
  393.  
  394. * pour le cas ou l'on a plusieurs noeuds
  395. * !!! hyp : Pnl de type POI1 !!!
  396. si (nPnl ega 1); geo1 = pnl;
  397. sinon; geo1 = manu 'SUPER' pnl;
  398. finsi;
  399. KNLFRE = MANU 'RIGIDITE' geo1 RIDEP LISTR1 'QUEL';
  400.  
  401. sinon;
  402. KNLFRE = VIDE 'RIGIDITE';
  403.  
  404. finsi;
  405.  
  406.  
  407. *
  408. *=== CALCUL DE CNL PAR LES LISTREELs ==================================
  409. *
  410. * cas ou il existe Cnl = -dfnl/d\dot x :
  411. si (FLVITE et (neg (dime CNLT) 0));
  412.  
  413. * preparation des objets resultats
  414. LISTR2 = prog;
  415. ntfr = (2**(OTFR-1)) + 1;
  416.  
  417. *----------------------- boucle sur les lignes de la Jacobienne CNL=-dfnl/dv
  418. * (= boucle sur les ddl temporels)
  419. REPE Bligne (DIME CNLT); iligne = &Bligne;
  420.  
  421. CNLT_i = CNLT . iligne;
  422. nCNLT_i = DIME CNLT_i;
  423.  
  424. * harmonique 0 : TFR(CNL=-dfnl/dv * 0) = 0 ==> CNLFR . * = 0
  425. CNLFR = TABLE; CNLFI = TABLE;
  426. * boucle sur les colonnes de la Jacobienne temporelle : dfnl_iligne/dv_icnl
  427. REPE Bcnl1 nCNLT_i; icnl = &Bcnl1;
  428. CNLFR . icnl = prog ntfr*0.;
  429. CNLFI . icnl = prog ntfr*0.;
  430. FIN Bcnl1;
  431.  
  432. * harmonique ik :
  433. * TFR(CNL * -kw sin(ik wt)) = w * (CNLCR.ik + i CNLCI.ik)
  434. * TFR(CNL * kw cos(ik wt)) = w * (CNLSR.ik + i CNLSI.ik)
  435. * tables de listreels
  436. CNLCR = TABLE; CNLCI = TABLE;
  437. CNLSR = TABLE; CNLSI = TABLE;
  438. ik = 0;
  439. REPE Bdtfr nhbm; ik = ik + 1;
  440. DFFCR = TABLE; DFFCI = TABLE;
  441. DFFSR = TABLE; DFFSI = TABLE;
  442. * boucle sur les colonnes de la Jacobienne temporelle : -dfnl_iligne/dv_icnl
  443. REPE Bcnl2 nCNLT_i; icnl = &Bcnl2;
  444. CNLTCik = (CNLT_i . icnl) * ( LISCOS2 . ik);
  445. CNLTSik = (CNLT_i . icnl) * ( LISSIN2 . ik);
  446. EVDFTC = EVOL MANU LISTDT CNLTCik;
  447. EVDFTS = EVOL MANU LISTDT CNLTSik;
  448. EVDFFC = TFR OTFR EVDFTC 'REIM';
  449. EVDFFS = TFR OTFR EVDFTS 'REIM';
  450. DFFCR . icnl = EXTR EVDFFC 'ORDO' 1;
  451. DFFCI . icnl = EXTR EVDFFC 'ORDO' 2;
  452. DFFSR . icnl = EXTR EVDFFS 'ORDO' 1;
  453. DFFSI . icnl = EXTR EVDFFS 'ORDO' 2;
  454. FIN Bcnl2;
  455. CNLCR . ik = DFFCR; CNLCI . ik = DFFCI;
  456. CNLSR . ik = DFFSR; CNLSI . ik = DFFSI;
  457. FIN Bdtfr;
  458.  
  459. * on range les valeurs dans LISTR2
  460. * -les n*(2H+1) premieres valeurs
  461. * --> ligne 1 : dF_0/dU = [ dF_0/dU_0 .. dF_0/dU_ik dF_0/dV_ik .. ]
  462. ik = 0;
  463. REPE BML1 nCNLT_i; icnl = &BML1;
  464. LISTR2 = LISTR2 ET ((EXTR CNLFR . icnl 1) );
  465. REPE BMLk nhbm; ik = &BMLk;
  466. LISTR2 = LISTR2 ET ((EXTR CNLCR . ik . icnl 1) );
  467. LISTR2 = LISTR2 ET ((EXTR CNLSR . ik . icnl 1) );
  468. FIN BMLk;
  469. * LISTR2 = LISTR2 ET ((EXTR CNLFR . icnl 1) * -1);
  470. * REPE BMLk nhbm; ik = &BMLk;
  471. * LISTR2 = LISTR2 ET ((EXTR CNLCR . ik . icnl 1) * -1);
  472. * LISTR2 = LISTR2 ET ((EXTR CNLSR . ik . icnl 1) * -1);
  473. * FIN BMLk;
  474. FIN BML1;
  475. *
  476. * -lignes (2,3) , (4,5) ... (cos,sin) ... : dF_ih/dU
  477. REPE BMKk nhbm; ih = &BMKk;
  478. * ligne cos : dF_ih/dU = [ dF_ih/dU_0 .. dF_ih/dU_ik dF_ih/dV_ik .. ]
  479. REPE BMKC1 nCNLT_i; icnl = &BMKC1;
  480. LISTR2 = LISTR2 ET ((EXTR CNLFR . icnl (ih + 1)) * 2 );
  481. REPE BMKCK nhbm; ik = &BMKCK;
  482. LISTR2 = LISTR2 ET
  483. ((EXTR CNLCR . ik . icnl (ih + 1)) * 2 );
  484. LISTR2 = LISTR2 ET
  485. ((EXTR CNLSR . ik . icnl (ih + 1)) * 2 );
  486. FIN BMKCK;
  487. * LISTR2 = LISTR2 ET ((EXTR CNLFR . icnl (ih + 1)) * -2 );
  488. * REPE BMKCK nhbm; ik = &BMKCK;
  489. * LISTR2 = LISTR2 ET
  490. * ((EXTR CNLCR . ik . icnl (ih + 1)) * -2 );
  491. * LISTR2 = LISTR2 ET
  492. * ((EXTR CNLSR . ik . icnl (ih + 1)) * -2 );
  493. * FIN BMKCK;
  494. FIN BMKC1;
  495. * ligne sin : dG_ih/dU = [ dG_ih/dU_0 .. dG_ih/dU_ik dG_ih/dV_ik .. ]
  496. REPE BMKS1 nCNLT_i; icnl = &BMKS1;
  497. LISTR2 = LISTR2 ET ((EXTR CNLFI . icnl (ih + 1)) * -2 );
  498. REPE BMKSK nhbm; ik = &BMKSK;
  499. LISTR2 = LISTR2 ET
  500. ((EXTR CNLCI . ik . icnl (ih + 1)) * -2 );
  501. LISTR2 = LISTR2 ET
  502. ((EXTR CNLSI . ik . icnl (ih + 1)) * -2 );
  503. FIN BMKSK;
  504. * LISTR2 = LISTR2 ET ((EXTR CNLFI . icnl (ih + 1)) * 2 );
  505. * REPE BMKSK nhbm; ik = &BMKSK;
  506. * LISTR2 = LISTR2 ET
  507. * ((EXTR CNLCI . ik . icnl (ih + 1)) * 2 );
  508. * LISTR2 = LISTR2 ET
  509. * ((EXTR CNLSI . ik . icnl (ih + 1)) * 2 );
  510. * FIN BMKSK;
  511. FIN BMKS1;
  512. FIN BMKk;
  513.  
  514.  
  515. FIN Bligne;
  516. *----------------------------------------- fin de boucle sur les lignes
  517.  
  518. si (non (exis RIDEP));
  519. RIDEP = MOTS;
  520. inl = 0;
  521. REPE Bcfnl1 (nCNLT_i / nPnl); inl = inl + 1;
  522. iCOMNL = EXTR COMP_FNL inl;
  523. RIDEP = RIDEP ET COMPU . iCOMNL ;
  524. FIN Bcfnl1;
  525. * pour le cas ou l'on a plusieurs noeuds
  526. * !!! hyp : Pnl de type POI1 !!!
  527. si (nPnl ega 1); geo1 = pnl;
  528. sinon; geo1 = manu 'SUPER' pnl;
  529. finsi;
  530. finsi;
  531. CNLFRE = MANU 'RIGIDITE' geo1 RIDEP (wrads * LISTR2) 'QUEL';
  532. * CNLFRE = MANU 'RIGIDITE' geo1 RIDEP (LISTR2) 'QUEL';
  533.  
  534. sinon;
  535. CNLFRE = VIDE 'RIGIDITE';
  536.  
  537. finsi;
  538.  
  539.  
  540. *=== par differentiation numerique ====================================
  541.  
  542. si IDIFF;
  543.  
  544. * on met localement a faux pour alleger l appel a charmeca
  545. IFKNL = FAUX;
  546. * petite valeur
  547. xpetit = 1.E-6 * (maxi DEPTOT 'ABS');
  548. si(xpetit < 1.E-10); xpetit=1.E-10; finsi;
  549. xm1 = 1. / xpetit;
  550. * rigidite de sortie
  551. KNLF1 = VIDE 'RIGIDITE';
  552.  
  553. *--------- boucle sur les ddls frequentiels =
  554. ncomp = dime RIDEP;
  555. * + boucle sur les noeuds
  556. repe bnoeud nPnl;
  557. si (ega (type pnl) 'POINT');
  558. pnl1 = pnl;
  559. sinon;
  560. pnl1 = poin pnl &bnoeud;
  561. finsi;
  562. * + boucle sur les composantes
  563. repe bcomp ncomp;
  564. comp1 = extr RIDEP &bcomp;
  565.  
  566. * perturbation du vecteur U
  567. * mess 'perturbe noeud #' (noeu pnl1) ' selon ' comp1 'de' xpetit;
  568. DUPP1 = MANU 'CHPO' pnl1 1 comp1 xpetit;
  569. DEPP1 = DEPTOT + DUPP1;
  570.  
  571. * passage temporel
  572. DEPTP1 VITTP1 = AFT1 DEPP1;
  573.  
  574. * fnl(u(t),v(t)) = ?
  575. TABAFT . 'DEP_NL' = DEPTP1;
  576. TABAFT . 'VIT_NL' = VITTP1;
  577. CHARMECA TABAFT;
  578. FNLTP1 = TABAFT . 'FNL_T';
  579.  
  580. * perturbation associee de FNL dans le domaine frequentiel
  581. FNLF1 = AFT2 FNLTP1;
  582. KNLF1 = KNLF1 et (-1.*xm1 *
  583. (MANU 'RIGI' (FNLF1 - FNLFRE) 'COLO' 'QUEL' pnl1 comp1) );
  584.  
  585. * * compariason des 2 approches :
  586. * mess '>>> differentiation:'; list (-1.*xm1 *FNLF1);
  587. * FNLF0 = (KNLFRE et CNLFRE) * DUPP1 * xm1;
  588. * mess '>>> TFR:'; list FNLF0;
  589.  
  590.  
  591. fin bcomp;
  592. fin bnoeud;
  593. *--------- fin de boucle sur les ddls frequentiels
  594.  
  595. IFKNL = VRAI;
  596.  
  597. finsi;
  598.  
  599.  
  600. FINSI;
  601. *=== FIN DU CALCUL DE KNL et CNL =======================================
  602.  
  603.  
  604.  
  605. ************************************************************************
  606. * STOCKAGE DE LA FORCE NON-LINEAIRE ET SA DERIVEE *
  607. ************************************************************************
  608.  
  609. TNL1 = TABL ;
  610. TNL1 . 'ADDI_SECOND' = FNLFRE;
  611. SI IFKNL;
  612. si IDIFF;
  613. TNL1 . 'ADDI_MATRICE' = KNLF1;
  614. sinon;
  615. TNL1 . 'ADDI_MATRICE' = KNLFRE et CNLFRE;
  616. finsi;
  617. FINSI;
  618.  
  619. * mess '--------------------------------------';
  620. * mess 'Fnl='; list FNLFRE;
  621. * * SI IFKNL;
  622. * mess '--------------------------------------';
  623. * mess 'KNL='; list KNLFRE;
  624. * mess '--------------------------------------';
  625. * * FINSI;
  626. * si (IFKNL et ((maxi FNLFRE abs) > 1.E-16)); erre 21; finsi;
  627.  
  628. FINPROC TNL1;
  629.  
  630.  
  631.  

© Cast3M 2003 - Tous droits réservés.
Mentions légales