Kalo te përmbajtja

Metoda Newton-Raphson

Linearizimi dhe rekurrenca iterative

Ekuacioni jolinear për zhvendosjen nyjore \(\boldsymbol{u}_{n+1}\) në kohën \(t_{n+1}\), i marrë te Puna virtuale e forcave të jashtme dhe montimi i ekuacionit global, zgjidhet me metodën Newton-Raphson. Zhvendosja nyjore \(\boldsymbol{u}_n\) deri në kohën \(t_n\) merret si e njohur dhe inkrementi i zhvendosjes \(\Delta\boldsymbol{u}\) si ndryshore e panjohur, për të përcaktuar

\[ \boldsymbol{u}_{n+1} = \boldsymbol{u}_n + \Delta\boldsymbol{u} \]

Më tej neglizhohet varësia e vektorit të forcës së jashtme nga zhvendosja nyjore dhe, duke marrë \(\boldsymbol{F}(\boldsymbol{u}_{n+1}) = \boldsymbol{F}_{n+1}\),

\[ \boldsymbol{Q}(\boldsymbol{u}_n + \Delta\boldsymbol{u}) - \boldsymbol{F}_{n+1} = \boldsymbol{0} \]

zgjidhet ekuacioni i mësipërm.

Në zgjidhjen aktuale \(\Delta\boldsymbol{u}\), përcaktohet ngurtësia tangjente

\[ \boldsymbol{K} = \left. \frac{\partial \boldsymbol{Q}}{\partial \boldsymbol{u}} \right|_{\boldsymbol{u}_n + \Delta\boldsymbol{u}} \]

Duke përdorur këtë për të linearizuar ekuacionin jolinear,

\[ \boldsymbol{K}\, d\boldsymbol{u} + \boldsymbol{Q}(\boldsymbol{u}_n + \Delta\boldsymbol{u}) - \boldsymbol{F}_{n+1} = \boldsymbol{0} \]

merret shprehja e mësipërme. Le të jetë \(d\boldsymbol{u}_i\) korrigjimi në iteracionin e \(i\)-të dhe vektori i mbetjes në fillim të iteracionit

\[ \boldsymbol{R}_{i-1} = \boldsymbol{F}_{n+1} - \boldsymbol{Q}(\boldsymbol{u}_n + \Delta\boldsymbol{u}) \]

atëherë rekurrenca iterative është

\[ \boldsymbol{K}_i\, d\boldsymbol{u}_i = \boldsymbol{R}_{i-1}, \qquad \Delta\boldsymbol{u} \leftarrow \Delta\boldsymbol{u} + d\boldsymbol{u}_i \]

në formën e mësipërme. Mbetja \(\boldsymbol{R}_i\) është madhësia që përfaqëson çekuilibrin e forcave ndaj gjendjes së ekuilibrit.

Ndërtimi i matricës tangjente të ngurtësisë

Ngurtësia tangjente \(\boldsymbol{K} = \partial\boldsymbol{Q}/\partial\boldsymbol{u}\) ndërtohet duke derivuar pjesërisht vektorin e forcës së brendshme të elementit të marrë te Diskretizimi i punës virtuale të forcave të brendshme ndaj zhvendosjes nyjore, duke integruar integrandët përkatës mbi domenin e çdo elementi dhe duke i montuar. Duke shënuar integrandin në nivel elementi me \(\boldsymbol{K}^e_X\) (notacioni i konfigurimit referencë, metoda TL) ose \(\boldsymbol{K}^e_x\) (notacioni i konfigurimit aktual, metoda UL), ngurtësia tangjente e elementit është

\[ \boldsymbol{K}^e = \int_{\Omega^e_0} \boldsymbol{K}^e_X\, dV \quad (\text{TL}), \qquad \boldsymbol{K}^e = \int_{\Omega^e} \boldsymbol{K}^e_x\, dv \quad (\text{UL}) \]

siç jepet më sipër. Më poshtë paraqiten format përfundimtare të integrandëve TL/UL. Në të dy rastet ato zbërthehen në shumën e një termi të ngurtësisë materiale (termi i zhvendosjes fillestare) dhe një termi të ngurtësisë gjeometrike (termi i sforcimit fillestar).

Metoda Total Lagrange

Në metodën Total Lagrange supozohet marrëdhënia lineare \(\dot{\boldsymbol{S}} = \boldsymbol{\mathsf{C}}:\dot{\boldsymbol{E}}\) ndërmjet shpejtësisë së sforcimit të dytë Piola-Kirchhoff \(\dot{\boldsymbol{S}}\) dhe shpejtësisë së deformimit Green-Lagrange \(\dot{\boldsymbol{E}}\). Kjo korrespondon me ligjet konstitutive të materialeve elastike lineare (materiale St. Venant-Kirchhoff) dhe materialeve hiperelastike; në FrontISTR metoda Total Lagrange përdoret për këto materiale. Në këtë rast, integrandi i ngurtësisë tangjente të elementit shkruhet në formë tensoriale si

\[ \delta\boldsymbol{u}^{eT}\, \boldsymbol{K}^e_X\, \dot{\boldsymbol{u}}^e = \dot{\boldsymbol{S}}:\delta\boldsymbol{E} + \boldsymbol{S}:(\delta\boldsymbol{F}^T \dot{\boldsymbol{F}}) \]

Termi i parë në anën e djathtë është termi i ngurtësisë materiale (termi i zhvendosjes fillestare), ndërsa termi i dytë është termi i ngurtësisë gjeometrike (termi i sforcimit fillestar).

Në zbatimin e FrontISTR, ky integrand llogaritet në formë matricore duke përdorur notacionin Voigt:

\[ \boldsymbol{K}^e_X = (\boldsymbol{B}_L + \boldsymbol{B}_{NL})^T\, \tilde{\boldsymbol{C}}\, (\boldsymbol{B}_L + \boldsymbol{B}_{NL}) + \boldsymbol{F}_9^T\, \boldsymbol{S}_9\, \boldsymbol{F}_9 \]

Matrica të veçanta janë si më poshtë. \(\boldsymbol{B}_L, \boldsymbol{B}_{NL}\) janë matricat B të paraqitura te Diskretizimi i punës virtuale të forcave të brendshme, ndërsa \(\tilde{\boldsymbol{C}}\) është matrica \(6\times 6\) e ngurtësisë materiale që përfaqëson tensorin konstitutiv \(\boldsymbol{\mathsf{C}}\) në notacionin Voigt (Notacioni tensorial dhe bazat matematikore). \(\boldsymbol{S}_9, \boldsymbol{F}_9\) janë matricat e rirenditjes të përdorura për të shprehur termin e ngurtësisë gjeometrike si prodhim matricor. Së pari, për një tensor të rendit të dytë \(3\times 3\) \(\boldsymbol{A}\), përkufizohet notacioni \([\,\cdot\,]\) që e rirendit në një vektor me 9 komponentë:

\[ [\boldsymbol{A}] = (A_{11}, A_{21}, A_{31}, A_{12}, A_{22}, A_{32}, A_{13}, A_{23}, A_{33})^T \]

Me këtë përkufizim, \(\boldsymbol{F}_9\) është matrica \(9\times d n_e\) që shpreh variacionin e gradientit të deformimit në formën \([\delta\boldsymbol{F}] = \boldsymbol{F}_9\, \delta\boldsymbol{u}^e\). Për nyjën e elementit \(\alpha = 1, \ldots, n_e\), blloku përkatës \(9\times d\) është

\[ [\boldsymbol{F}_9]_\alpha = \begin{bmatrix} (\partial N_\alpha^e/\partial X_1)\, \boldsymbol{I} \\ (\partial N_\alpha^e/\partial X_2)\, \boldsymbol{I} \\ (\partial N_\alpha^e/\partial X_3)\, \boldsymbol{I} \end{bmatrix} \qquad (\boldsymbol{I} \text{ është matrica njësi } 3\times 3) \]

dhe \(\boldsymbol{F}_9 = [[\boldsymbol{F}_9]_1, \ldots, [\boldsymbol{F}_9]_{n_e}]\) merret duke i vendosur këto blloqe horizontalisht sipas rendit të nyjeve të elementit. \(\boldsymbol{S}_9\) është matrica \(9\times 9\) e zgjedhur në mënyrë që, kur kombinohet me këtë matricë, termi i ngurtësisë gjeometrike të shprehet si \(\delta\boldsymbol{u}^{eT}\, \boldsymbol{F}_9^T \boldsymbol{S}_9 \boldsymbol{F}_9\, \dot{\boldsymbol{u}}^e\):

\[ \boldsymbol{S}_9 = \begin{bmatrix} S_{11} \boldsymbol{I} & S_{12} \boldsymbol{I} & S_{13} \boldsymbol{I} \\ S_{21} \boldsymbol{I} & S_{22} \boldsymbol{I} & S_{23} \boldsymbol{I} \\ S_{31} \boldsymbol{I} & S_{32} \boldsymbol{I} & S_{33} \boldsymbol{I} \end{bmatrix} \]

Kjo është matrica përkatëse.

Metoda Updated Lagrange

Në metodën Updated Lagrange supozohet marrëdhënia lineare \(\hat{\boldsymbol{\sigma}}^{\nabla J} = \boldsymbol{\mathsf{C}}:\boldsymbol{D}\) ndërmjet shpejtësisë Jaumann të tensorit relativ të sforcimit Kirchhoff \(\hat{\boldsymbol{\sigma}}^{\nabla J}\) dhe tensorit të shpejtësisë së deformimit \(\boldsymbol{D}\). Kjo është forma e ligjit konstitutiv hipoelastik të përbashkët për materialet elastike lineare, elastoplastike dhe me zvarritje; në FrontISTR metoda Updated Lagrange përdoret për këto materiale. Në këtë rast, integrandi i ngurtësisë tangjente të elementit, i shprehur në konfigurimin aktual, shkruhet në formë tensoriale si

\[ \delta\boldsymbol{u}^{eT}\, \boldsymbol{K}^e_x\, \dot{\boldsymbol{u}}^e = \boldsymbol{\sigma}^{\nabla T}:\delta\boldsymbol{A}_{(L)} + \boldsymbol{\sigma}:(\delta\boldsymbol{F}_t^T\, \boldsymbol{L}) \]

ku \(\boldsymbol{\sigma}^{\nabla T}\) është shpejtësia Truesdell, \(\boldsymbol{A}_{(L)}\) pjesa lineare e deformimit Almansi, \(\boldsymbol{F}_t = \partial\boldsymbol{u}/\partial\boldsymbol{x}\) gradienti i zhvendosjes ndaj konfigurimit aktual dhe \(\boldsymbol{L}\) tensori i gradientit të shpejtësisë. Termi i parë në anën e djathtë është termi i ngurtësisë materiale, ndërsa i dyti është termi i ngurtësisë gjeometrike.

Në zbatimin e FrontISTR, ky integrand llogaritet në formë matricore duke përdorur notacionin Voigt:

\[ \boldsymbol{K}^e_x = \boldsymbol{b}^T\, (\tilde{\boldsymbol{C}} - \boldsymbol{G})\, \boldsymbol{b} + \boldsymbol{f}_9^T\, \boldsymbol{\sigma}_9\, \boldsymbol{f}_9 \]

\(\boldsymbol{b}\) është matrica B e ndërtuar në konfigurimin aktual (Diskretizimi i punës virtuale të forcave të brendshme). \(\boldsymbol{\sigma}_9, \boldsymbol{f}_9\) merren nga \(\boldsymbol{S}_9, \boldsymbol{F}_9\) të përkufizuara për metodën TL duke zëvendësuar sforcimin e dytë PK \(\boldsymbol{S}\) me sforcimin Cauchy \(\boldsymbol{\sigma}\) dhe gradientin e konfigurimit referencë \(\partial N_\alpha^e/\partial X_i\) me gradientin e konfigurimit aktual \(\partial N_\alpha^e/\partial x_i\).

\(\boldsymbol{G}\) është një matricë korrigjimi e varur nga sforcimi Cauchy, e nevojshme për ta bërë ligjin konstitutiv hipoelastik \(\hat{\boldsymbol{\sigma}}^{\nabla J} = \boldsymbol{\mathsf{C}}:\boldsymbol{D}\) të përputhshëm me kuadrin e ngurtësisë tangjente si ligj konstitutiv i bazuar në shpejtësinë Truesdell. Ajo merret duke vendosur komponentët e tensorit të rendit të katërt \(G_{ijkl} = \delta_{il}\sigma_{kj} + \delta_{jl}\sigma_{ik}\) në notacionin Voigt \(6\times 6\):

\[ \boldsymbol{G} = \begin{bmatrix} 2\sigma_{11} & 0 & 0 & \sigma_{12} & 0 & \sigma_{31} \\ 0 & 2\sigma_{22} & 0 & \sigma_{12} & \sigma_{23} & 0 \\ 0 & 0 & 2\sigma_{33} & 0 & \sigma_{23} & \sigma_{31} \\ \sigma_{12} & \sigma_{12} & 0 & \tfrac{\sigma_{11}+\sigma_{22}}{2} & \tfrac{\sigma_{12}}{2} & \tfrac{\sigma_{23}}{2} \\ 0 & \sigma_{23} & \sigma_{23} & \tfrac{\sigma_{12}}{2} & \tfrac{\sigma_{22}+\sigma_{33}}{2} & \tfrac{\sigma_{12}}{2} \\ \sigma_{31} & 0 & \sigma_{31} & \tfrac{\sigma_{23}}{2} & \tfrac{\sigma_{12}}{2} & \tfrac{\sigma_{33}+\sigma_{11}}{2} \end{bmatrix} \]

Kjo është matrica përkatëse.

Montimi i matricës globale të ngurtësisë

Ngurtësia tangjente globale \(\boldsymbol{K}\) merret duke ndarë ngurtësinë e elementit \(\boldsymbol{K}^e\) në blloqe \(d\times d\) \(\boldsymbol{K}^e_{\alpha\beta}\) për çdo çift nyjesh dhe duke përdorur bashkësinë e montimit për tensorë të rendit të dytë \(\mathcal{E}^2(i_g, i_h)\), të paraqitur te Montimi i madhësive fizike nyjore të elementeve:

\[ \boldsymbol{K}_{i_gi_h} = \sum_{(e,\alpha,\beta) \in \mathcal{E}^2(i_g, i_h)} \boldsymbol{K}^e_{\alpha\beta} \]

Vlerat vendosen në një matricë me rreshtin \(i_g\) dhe kolonën \(i_h\). Në zbatim, bashkësia \(\mathcal{E}^2\) nuk ndërtohet shprehimisht; blloqet përkatëse shtohen drejtpërdrejt brenda ciklit të elementeve. Matrica është katrore me dimensionin shkallë lirie për nyje \(\times\) numri i përgjithshëm i nyjeve \(n_g\), por meqenëse komponentët ndërmjet nyjeve që nuk lidhen përmes elementeve janë \(0\), ajo ruhet në format matrice të rrallë.

Matricat e ngurtësisë së elementit për metodat TL dhe UL kanë të njëjtën formë, përveç ndryshimit të konfigurimit referencë (koordinatat nyjore dhe burimi për ndërtimin e matricës B) dhe pranisë ose mungesës së matricës \(\boldsymbol{G}\). Prandaj, FrontISTR i zbaton të dyja me një nënprogram të përbashkët.

Algoritmi i iteracionit

Duke përmbledhur sa më sipër, në fillim të iteracionit vendoset \(\Delta\boldsymbol{u} = \boldsymbol{0}\) dhe llogaritet mbetja fillestare \(\boldsymbol{R}_0 = \boldsymbol{F}_{n+1} - \boldsymbol{Q}(\boldsymbol{u}_n)\). Pastaj, në iteracionin e \(i\)-të kryhen hapat vijues.

  1. Në zhvendosjen aktuale \(\boldsymbol{u}_n + \Delta\boldsymbol{u}\), llogaritni ngurtësinë tangjente \(\boldsymbol{K}_i = \partial\boldsymbol{Q}/\partial\boldsymbol{u}|_{\boldsymbol{u}_n + \Delta\boldsymbol{u}}\) sipas procedurës te Ndërtimi i matricës tangjente të ngurtësisë.
  2. Për të zbatuar kushtet kufitare gjeometrike, modifikoni matricën tangjente të ngurtësisë dhe vektorin e mbetjes për shkallët e lirisë ku janë vendosur kufizime zhvendosjeje, duke marrë \(\tilde{\boldsymbol{K}}_i, \tilde{\boldsymbol{R}}_{i-1}\) (Trajtimi i kushteve kufitare gjeometrike).
  3. Zgjidhni ekuacionin linear \(\tilde{\boldsymbol{K}}_i\, d\boldsymbol{u}_i = \tilde{\boldsymbol{R}}_{i-1}\) për të marrë korrigjimin \(d\boldsymbol{u}_i\). Kjo procedurë shpesh zë pjesën më të madhe të kostos llogaritëse të iteracionit.
  4. Përditësoni inkrementin e zhvendosjes si \(\Delta\boldsymbol{u} \leftarrow \Delta\boldsymbol{u} + d\boldsymbol{u}_i\) dhe, në përputhje me këtë, llogaritni vektorin e forcës së brendshme \(\boldsymbol{Q}(\boldsymbol{u}_n + \Delta\boldsymbol{u})\) dhe mbetjen \(\boldsymbol{R}_i = \boldsymbol{F}_{n+1} - \boldsymbol{Q}(\boldsymbol{u}_n + \Delta\boldsymbol{u})\).
  5. Kontrolloni konvergjencën dhe përfundoni iteracionin nëse ajo është arritur. Në komponentët e mbetjes \(\boldsymbol{R}_i\) që u përkasin shkallëve të lirisë me kushte kufitare gjeometrike shfaqen komponentë që korrespondojnë me reaksionet e kufizimit; prandaj treguesi i konvergjencës ndërtohet nga \(\tilde{\boldsymbol{R}}_i\) pasi këta komponentë përjashtohen. Treguesit dhe pragjet konkrete të konvergjencës trajtohen te Kriteret e konvergjencës. Nëse nuk arrihet konvergjenca dhe arrihet kufiri i iteracioneve, iteracioni konsiderohet i dështuar.

Kur iteracioni konvergon, \(\Delta\boldsymbol{u}\) i konverguar i shtohet \(\boldsymbol{u}_n\) për të marrë zhvendosjen e akumuluar në kohën \(t_{n+1}\), \(\boldsymbol{u}_{n+1} = \boldsymbol{u}_n + \Delta\boldsymbol{u}\), dhe vazhdohet me hapin kohor pasues.

Tema të lidhura

AI-assisted translation May contain errors Official docs Status