Агуулгыг алгасах

Ньютон-Рафсоны арга

Шугамчлах ба давталтын рекуррент хамаарал

Гадаад хүчний виртуал ажил ба глобал тэгшитгэлийн угсралт-аас \(t_{n+1}\) хугацаан дахь зангилааны шилжилт \(\boldsymbol{u}_{n+1}\)-ийн шугаман бус тэгшитгэл гарч, үүнийг Ньютон-Рафсоны аргаар шийднэ. \(t_n\) хугацаа хүртэлх зангилааны шилжилт \(\boldsymbol{u}_n\) мэдэгдсэн гэж үзээд, тодорхойлох үл мэдэгдэх хувьсагчаар шилжилтийн өсөлт \(\Delta\boldsymbol{u}\)-г авна

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

Цаашид гадаад хүчний векторын зангилааны шилжилтээс хамаарах хамаарлыг үл тооцож, \(\boldsymbol{F}(\boldsymbol{u}_{n+1}) = \boldsymbol{F}_{n+1}\) гэж авбал

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

тэгшитгэлийг шийднэ.

Одоогийн шийдэл \(\Delta\boldsymbol{u}\) дээр шүргэгч хөшүүн чанар-ыг

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

Үүнийг ашиглан шугаман бус тэгшитгэлийг шугамчилбал

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

\(i\)-р давталтын засварыг \(d\boldsymbol{u}_i\), давталт эхлэх үеийн үлдэгдлийн векторыг

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

Тэгвэл давталтын рекуррент хамаарал

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

Иймээс үлдэгдэл \(\boldsymbol{R}_i\) нь тэнцвэрээс зөрсөн хүчний тэнцвэргүй байдлыг илэрхийлэх хэмжигдэхүүн юм.

Шүргэгч хөшүүн чанарын матрицыг байгуулах нь

Шүргэгч хөшүүн чанар \(\boldsymbol{K} = \partial\boldsymbol{Q}/\partial\boldsymbol{u}\)-ыг Дотоод хүчний виртуал ажлын дискретчлэл-д гаргасан элементийн дотоод хүчний векторыг зангилааны шилжилтийн хувьд хэсэгчлэн уламжилж, гарсан элементийн түвшний интегралын илэрхийллийг элемент бүрийн муж дээр интегралчлан угсрах замаар байгуулна. Элементийн түвшний интегралын илэрхийллийг \(\boldsymbol{K}^e_X\) (эхний байрлалын тэмдэглэгээ, TL томьёолол) эсвэл \(\boldsymbol{K}^e_x\) (одоогийн байрлалын тэмдэглэгээ, UL томьёолол) гэж тэмдэглэвэл элементийн шүргэгч хөшүүн чанар

\[ \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}) \]

Доор TL/UL интегралын илэрхийллүүдийн эцсийн хэлбэрийг үзүүлэв. Аль алинд нь материалын хөшүүн чанарын гишүүн (анхны шилжилтийн гишүүн) болон геометрийн хөшүүн чанарын гишүүн (анхны хүчдэлийн гишүүн)-ийн нийлбэр болгон задална.

Total Lagrange томьёолол

Total Lagrange томьёололд хоёрдугаар Piola-Kirchhoff хүчдэлийн хурд \(\dot{\boldsymbol{S}}\) ба Green-Lagrange харьцангуй деформацийн хурд \(\dot{\boldsymbol{E}}\)-ийн хооронд \(\dot{\boldsymbol{S}} = \boldsymbol{\mathsf{C}}:\dot{\boldsymbol{E}}\) шугаман хамаарал байна гэж үзнэ. Энэ нь шугаман уян материал (St. Venant-Kirchhoff материал) болон хэт уян харимхай материалын бүтцийн хуулиудад тохирох бөгөөд FrontISTR эдгээр материалд Total Lagrange томьёоллыг хэрэглэнэ. Тэгвэл элементийн шүргэгч хөшүүн чанарын интегралын илэрхийллийг тензор хэлбэрээр

\[ \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}}) \]

Баруун талын эхний гишүүн нь материалын хөшүүн чанарын гишүүн (анхны шилжилтийн гишүүн), хоёр дахь гишүүн нь геометрийн хөшүүн чанарын гишүүн (анхны хүчдэлийн гишүүн) юм.

FrontISTR-ийн хэрэгжүүлэлтэд энэ интегралын илэрхийллийг 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 \]

Матриц бүрийг дараах байдлаар тодорхойлно. \(\boldsymbol{B}_L, \boldsymbol{B}_{NL}\) нь Дотоод хүчний виртуал ажлын дискретчлэл-д танилцуулсан B матрицууд, \(\tilde{\boldsymbol{C}}\) нь бүтцийн тензор \(\boldsymbol{\mathsf{C}}\)-ийн Voigt дүрслэл буюу \(6\times 6\) материалын хөшүүн чанарын матриц (Тензор тэмдэглэгээ ба математик үндэс) юм. \(\boldsymbol{S}_9, \boldsymbol{F}_9\) нь геометрийн хөшүүн чанарын гишүүнийг матрицын үржвэрээр илэрхийлэхэд хэрэглэгдэх дараах дахин байрлуулах матрицууд. Эхлээд \(3\times 3\) хоёрдугаар эрэмбийн тензор \(\boldsymbol{A}\)-г 9 бүрэлдэхүүнт вектор болгон дахин байрлуулах \([\,\cdot\,]\) тэмдэглэгээг

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

Энэ тодорхойлолтоор \(\boldsymbol{F}_9\) нь деформацийн градиентийн вариацыг \([\delta\boldsymbol{F}] = \boldsymbol{F}_9\, \delta\boldsymbol{u}^e\) хэлбэрээр илэрхийлэх \(9\times d n_e\) матриц байна. Элементийн зангилаа \(\alpha = 1, \ldots, n_e\) бүрт харгалзах \(9\times d\) блок нь

\[ [\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{ нь } 3\times 3 \text{ нэгж матриц}) \]

ба блокуудыг элементийн зангилааны дарааллаар хэвтээ байрлуулсан \(\boldsymbol{F}_9 = [[\boldsymbol{F}_9]_1, \ldots, [\boldsymbol{F}_9]_{n_e}]\) гэж өгнө. \(\boldsymbol{S}_9\)-г энэ матрицтай хослуулахад геометрийн хөшүүн чанарын гишүүн \(\delta\boldsymbol{u}^{eT}\, \boldsymbol{F}_9^T \boldsymbol{S}_9 \boldsymbol{F}_9\, \dot{\boldsymbol{u}}^e\) хэлбэрээр илэрхийлэгдэхээр сонгох бөгөөд дараах \(9\times 9\) матриц байна

\[ \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} \]

Ингэж матриц гарна.

Updated Lagrange томьёолол

Updated Lagrange томьёололд харьцангуй Кирхгофын хүчдэлийн тензорын Жауманны хурд \(\hat{\boldsymbol{\sigma}}^{\nabla J}\) ба деформацийн хурдны тензор \(\boldsymbol{D}\)-ийн хооронд \(\hat{\boldsymbol{\sigma}}^{\nabla J} = \boldsymbol{\mathsf{C}}:\boldsymbol{D}\) шугаман хамаарал байна гэж үзнэ. Энэ нь шугаман уян, уян-хуванцар болон мөлхөлтийн материалд нийтлэг гипоуян бүтцийн хуулийн хэлбэр бөгөөд FrontISTR эдгээр материалд Updated Lagrange томьёоллыг хэрэглэнэ. Одоогийн байрлалд илэрхийлсэн элементийн шүргэгч хөшүүн чанарын интегралын илэрхийллийг тензор хэлбэрээр

\[ \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}) \]

энд \(\boldsymbol{\sigma}^{\nabla T}\) нь Truesdell хурд, \(\boldsymbol{A}_{(L)}\) нь Almansi харьцангуй деформацийн шугаман хэсэг, \(\boldsymbol{F}_t = \partial\boldsymbol{u}/\partial\boldsymbol{x}\) нь одоогийн байрлалын хувьд авсан шилжилтийн градиент, \(\boldsymbol{L}\) нь хурдны градиентийн тензор. Баруун талын эхний гишүүн нь материалын хөшүүн чанарын гишүүн, хоёр дахь нь геометрийн хөшүүн чанарын гишүүн.

FrontISTR-ийн хэрэгжүүлэлтэд энэ интегралын илэрхийллийг 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}\) нь одоогийн байрлалд байгуулагдсан B матриц (Дотоод хүчний виртуал ажлын дискретчлэл). \(\boldsymbol{\sigma}_9, \boldsymbol{f}_9\)-г TL томьёололд тодорхойлсон \(\boldsymbol{S}_9, \boldsymbol{F}_9\)-д хоёрдугаар PK хүчдэл \(\boldsymbol{S}\)-ийг Cauchy хүчдэл \(\boldsymbol{\sigma}\)-аар, эхний байрлалын градиент \(\partial N_\alpha^e/\partial X_i\)-г одоогийн байрлалын градиент \(\partial N_\alpha^e/\partial x_i\)-ээр сольж олно.

\(\boldsymbol{G}\) нь гипоуян бүтцийн хууль \(\hat{\boldsymbol{\sigma}}^{\nabla J} = \boldsymbol{\mathsf{C}}:\boldsymbol{D}\)-г Truesdell хурд дээр суурилсан бүтцийн хууль хэлбэрийн шүргэгч хөшүүн чанарын хүрээтэй нийцүүлэхэд шаардлагатай, Cauchy хүчдэлээс хамаарах засварын матриц юм. Үүнийг дөрөвдүгээр эрэмбийн тензорын бүрэлдэхүүн \(G_{ijkl} = \delta_{il}\sigma_{kj} + \delta_{jl}\sigma_{ik}\)\(6\times 6\) Voigt хэлбэрээр байрлуулж

\[ \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} \]

Ингэж матриц гарна.

Глобал хөшүүн чанарын матрицыг угсрах нь

Глобал шүргэгч хөшүүн чанар \(\boldsymbol{K}\)-г элементийн хөшүүн чанар \(\boldsymbol{K}^e\) бүрийг зангилааны хос бүрийн \(d\times d\) блок \(\boldsymbol{K}^e_{\alpha\beta}\) болгон хувааж, Элемент-зангилааны физик хэмжигдэхүүний угсралт-д танилцуулсан хоёрдугаар эрэмбийн тензорын угсралтын олонлог \(\mathcal{E}^2(i_g, i_h)\)-г ашиглан олно:

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

Гарсан утгуудыг \(i_g\) мөр, \(i_h\) баганатай матриц болгон байрлуулна. Хэрэгжүүлэлтэд \(\mathcal{E}^2\) олонлогийг илээр байгуулдаггүй; харин харгалзах блокуудыг элементийн давталт дотор шууд нэмнэ. Матриц нь нэг зангилааны чөлөөний зэрэг \(\times\) нийт зангилааны тоо \(n_g\) хэмжээтэй квадрат матриц боловч элементээр холбогдоогүй зангилаануудын хоорондох бүрэлдэхүүн \(0\) тул сийрэг матриц хэлбэрээр хадгална.

TL ба UL томьёоллын элементийн хөшүүн чанарын матрицууд нь эхний байрлалыг солих (зангилааны координат ба B матриц байгуулах эх үүсвэр) болон \(\boldsymbol{G}\) матриц байгаа эсэхээс бусдаар ижил хэлбэртэй. Иймээс FrontISTR хоёр томьёоллыг нэг нийтлэг дэд программд хэрэгжүүлдэг.

Давталтын алгоритм

Дээрхийг нэгтгэн, давталт эхлэхэд \(\Delta\boldsymbol{u} = \boldsymbol{0}\) гэж тавьж анхны үлдэгдэл \(\boldsymbol{R}_0 = \boldsymbol{F}_{n+1} - \boldsymbol{Q}(\boldsymbol{u}_n)\)-г тооцоолно. Дараа нь \(i\)-р давталтад дараах үйлдлүүдийг гүйцэтгэнэ.

  1. Одоогийн шилжилт \(\boldsymbol{u}_n + \Delta\boldsymbol{u}\) дээр Шүргэгч хөшүүн чанарын матрицыг байгуулах хэсгийн аргаар шүргэгч хөшүүн чанар \(\boldsymbol{K}_i = \partial\boldsymbol{Q}/\partial\boldsymbol{u}|_{\boldsymbol{u}_n + \Delta\boldsymbol{u}}\)-г тооцоолно.
  2. Геометрийн хилийн нөхцөлийг хэрэгжүүлэхийн тулд шилжилтийн хязгаарлалттай чөлөөний зэргийн шүргэгч хөшүүн чанарын матриц ба үлдэгдлийн векторыг өөрчилж \(\tilde{\boldsymbol{K}}_i, \tilde{\boldsymbol{R}}_{i-1}\)-г олно (Геометрийн хилийн нөхцөлийг боловсруулах).
  3. Шугаман тэгшитгэл \(\tilde{\boldsymbol{K}}_i\, d\boldsymbol{u}_i = \tilde{\boldsymbol{R}}_{i-1}\)-г шийдэж засвар \(d\boldsymbol{u}_i\)-г олно. Энэ үйлдэл нь давталтын тооцооллын зардлын ихэнхийг эзлэх нь элбэг.
  4. Шилжилтийн өсөлтийг \(\Delta\boldsymbol{u} \leftarrow \Delta\boldsymbol{u} + d\boldsymbol{u}_i\) гэж шинэчилж, түүнд харгалзан дотоод хүчний вектор \(\boldsymbol{Q}(\boldsymbol{u}_n + \Delta\boldsymbol{u})\) болон үлдэгдэл \(\boldsymbol{R}_i = \boldsymbol{F}_{n+1} - \boldsymbol{Q}(\boldsymbol{u}_n + \Delta\boldsymbol{u})\)-г тооцоолно.
  5. Нийлэлтийг шалгаж, нийлэлсэн бол давталтыг дуусгана. Геометрийн хилийн нөхцөлтэй чөлөөний зэрэгт үлдэгдэл \(\boldsymbol{R}_i\) дотор хязгаарлагдсан урвалд харгалзах бүрэлдэхүүнүүд гардаг тул эдгээрийг хассан \(\tilde{\boldsymbol{R}}_i\)-аас нийлэлтийн үзүүлэлтийг байгуулна. Тодорхой нийлэлтийн үзүүлэлт ба босгыг Нийлэлтийн шалгуур-т тайлбарласан. Нийлэлтэд хүрээгүй бөгөөд давталтын хязгаарт хүрсэн бол давталтыг амжилтгүй гэж үзнэ.

Давталт нийлэлсэн үед нийлэлсэн \(\Delta\boldsymbol{u}\)\(\boldsymbol{u}_n\)-д нэмж \(t_{n+1}\) хугацаан дахь хуримтлагдсан шилжилт \(\boldsymbol{u}_{n+1} = \boldsymbol{u}_n + \Delta\boldsymbol{u}\)-г олно, дараагийн хугацааны алхамд шилжинэ.

Холбогдох сэдвүүд