Смагоринскийн загварчлал /Smagorinsky Sub-grid scale model/

Байгаль дээрх ихэнх урсгал турбулент шинж чанартай байдаг. Голын усны урсгал, шуурга, агаарын урсгал, циклон, үүлний хэв маягийн өөрчлөлт эсвэл нарны шуурга, огторгуй дахь урсгал гэх мэт гээд турбулент шинж чанарыг агуулсан урсгал эргэн тойронд элбэг бий. Өмнө бид К-Эпсилон загварын тухай авч үзсэн бол одоо Том хуйлралтын симуляцийн нэг хувилбар болох Смагорискийн дэд торны ойролцооллыг авч үзэх гэж байна. Дифференциал тэгшитгэлийг тооцон бодоход бодлогын хүрээг тодорхой геометрийн элементүүдэд торлон хуваах ба элементийн уртын хэмжээс нь тооцооны чухал парамет болдог. Уг уртын хэмжээс дээр том хуйлралтыг ялган таних дахин нэг уртын хэмжээсийг авч үзэж байгаа нь дэд тор гэж нэрлэгдэж байгаа юм. 

Даралт-Хурдны итерацитай Шийдлийн алгоритмоор НС-ын тэгшитгэлийг боловсруулсан код. Турбулент Смагорискийн загвар мөн ашиглагдсан. Нягтын үйлчлэлийн бодлого /Lock-exchange by Prof Hosoyamada Tokuzo/

Смагоринскийн дэд торны хэмжээст загварчлал

Смагоринскийн турбулент загварчлал нь Том хуйлралтын симуляци (Large Eddy Simulation) хэмээх турбулент урсгалыг загварчлах багц загварын нэг хэлтэрхий юм. 1963 онд Жосеф Смагоринский атмосферийн эрчимжсэн урсгалыг тооцолохдоо том хуйлралтын симуляцийг танилцуулсан боловч анх 1970 онд Жамес Леардорф үндсэн санааг нь боловсруулсан байна.

Шингэний хөдөлгөөн нь орон зай болон цаг хугацааны хэмжээстэд турбулент урсгалыг тооцоолон илэрхийлэх боломжтой ба энд зарах аргачлалыг Шууд тооцооллын арга (directnumerical simulation) гэж нэрлэдэг. Энэ аргачлал нь тооцооллын багтаамж өндөртэй, урт хугацааны симуляци шаарддаг учир цаг хугацаа болон багтаамж хожих турбулент загваруудыг боловсруулж эхэлсэн байна.

Нөгөө талаас Колмогоровийн (1941) хувийн төсийн (self-similarity) онолын утга санаагаар бол урсгалын том хуй нь бага хэмжээсийн хувьд нийтлэг мэт болж геометрын хэмжээсээс хамааралтай болно.  Энэ онцлог санаан дээр Том хуйлралтын симуляцийн загварууд боловсруулагдсан байна. Ингээд Смагоринскийн дэд торны хэмжээст загварчлалын тухай авч үзье.

Шингэний динамикийг дүрслэн харуулах тэгшитгэл нь урсгал тасралтгүйн болон моментын тэгшитгэлүүд билээ. Тэдгээрийг үл шахагдах шингэний хувьд бичвэл:

Энд тау’’ нь зунгааралтын шүргэх хүчдэлийн тенсор ба f_m нь биеийн хүч болно.

Тэгшитгэлийг шүүх /фильтр/

Том хуйлралтын тооцооллын гол зорилго нь том хэмжээст турбулент урсгалыг жижиг хэмжээст загварчлал ашиглан ил байдлаар тооцох явдал юм. Үүнийг хамгийн сайн гүйцэтгэхэд шүүгдсэн тэгшитгэл хэрэгтэй болно. Нэг хэмжээст тэмдэглэгээг ашиглан хурдны шүүлтийг илэрхийлбэл:

Энд G(x,x’) нь шүүлтүүр функц бөгөөд хэд хэдэн хэлбэртэй байх орчны функц юм. Тухайлбал Гауссын шүүлтүүр, Хайрцаган шүүлтүүр, таслах шүүлтүүр гэх мэт олон шүүлтүүрийн цөм функц бий. Шүүлтүүр бүр уртын хэмжээс гэж нэрлэгдэх gulwaljin параметертэй холбогдоно. Ерөнхий тохиолдолд уртын хэмжээсээс илүү гарах хуйлралтууд нь том хуйлралтын ангилалд багтаж дахин бодогдох хэмжээст хамаарагдах ба түүнээс бага хэмжээтэй нь жижиг хуйлралтанд багтаж тооцогдсон байх хэрэгтэй. Шүүлтүүрийг хэрэглэсэний дараа хурдыг илэрхийлбэл:

Үүнд u_i нь дахин шийдэгдэх хэмжээс дахь хэсэг, u_i’ нь дэд торны хэмжээст хэсэг болно. Ерөнхийдөө практик талдаа аргачлал нь хэлхээсийг өөрийг (хайрцаган шүүлтүүр) нь шүүлтүүр болгон ашиглах явдал юм. Дээрхтэй адил даралтын хувьд шүүлтүүр болон хэсэгчилсэн даралтын илэрхийллийг бичиж болно.

НС-ын тэгшитгэлд шүүлтүүрийг оруулах үед гарах тэгшитгэлийн багц нь Рейнольдсын Дундаж НС-ын тэгшитгэлтэй ижил хэлбэртэй байна.

Тэгшитгэл 5-аас ажиглавал моментийн тэгшитгэлийн шугаман бус байдал нь РДНС-ын тэгшитгэл дахь Рейнольдсын хүчдэлтэй адил гишүүнийг бий болгох ба үүнийг РДНС тэгшитгэл талаас авч үзвэл:

байх бөгөөд тэнцэтгэл бусын зүүн гар талын илэрхийлэл нь тооцоологдох боломжгүй. Энэ тэнцэтгэл бусын балансаас дараах илэрхийлэл гарч ирнэ.

Энд Tau^s_ij нь дэд хэмжээс дахь (SGS) Рейнольдсын хүчдэл юм. Үүнийг загварчлах ёстой. Энэ бол физик талаасаа хүчдэл биш бөгөөд жижиг эсвэл дахин шийдэгдээгүй хэмжээсийн улмаар бий болох том хэмжээст моментын цүлхэлт гэж үзэж болно. Тэгшитгэл 8-г тэгшитгэл 6-руу орлуулж дараагийн тэгшитгэлийг гарган авбал:

Тэгшитгэл 9 дахь tau^s_ij-ийн ойролцоолол нь Том хуйлралтын симуляцийн гол асуудал байх ба энэ нь өөр өөр төрлийн том хуйлралтын загварчлалын эх сурвалж болдог. Жишээ нь Смагорински-Лилигийн загвар, Алгебр-динамикийн загвар, Динамик-нийтлэг коэффицентын загвар, орчны динамик загвар, Дахин нормчилсон том хуйлралтын загвар, байгууламжийн загварчлал гэх мэт байна. Энд бид Смагоринскийн загварчлалыг авч үзье. Үүнийг дараах хэсгээр тайлбарлая.

Смагоринскийн дөхөлт

Энэ загварыг 1963 онд Смагорински санал болгосон бөгөөд өнөө үеийн хамгийн түгээмэл хэрэглэгдэх загварын нэг болж байна. Зарим судлаачид үүнийг Прандтлийн холигдох замын уртын загварын том хуйлралтын симуляци дахь хувилбар хэмээн үздэг. Ламинар урсгал дахь хүчдэлийн нөлөөтэй ижил болгож Дэд хэмжээс дахь /дэд торны/ хүчдэлийг бичвэл:

Энд ню_Т нь хуйлралтын зунгааралт, S^bar_ij нь дахин тооцогдсон хурдны орон дахь сунах хүчдэлийн утгыг илэхрхийлэх ба дараах байдлаар тодорхойлогдоно.

Хэмжээсийн шинжилгээг хэрэглэж анализ хийвэл хуйлралтын зунгааралтын тохиромжит хэлбэр нь дараах байдалтай байх нь батлагдана.

Параметр С_S нь изотрофик турбулентлэг чанарыг агуулах өөр өөр онолоор тодорхойлогдох боломжтой ба тухайлбал 0.1-0.2-той ойролцоогоор авч үзэж болно. Гэвч энэ параметр нь Рейнольдсын тоо гэх мэт бусад хүчин зүйлүүдээс хамаарах функц гэдгийг мартаж болохгүй. Тухайлбал суваг дахь шингэний урсгалыг тооцоход энэ параметр нь 0.2-0.065 хүртэл буурч тодорхойлогдох тохиолдол бий. Энэ зөрүүтэй утгыг ашиглахад хуйлралтын зунгааралтын утга нэг эрэмбээр буурдаг байна. Нэмж хэлэхэд хананд ойрхон бүсэд параметрийн утгыг өшөө багасгаж өгөх шаардлага гардаг. РДНС-ын загвар дахь ханатай ойролцоо хуйлралтын зунгааралтыг багасгахад хамгийн тохиромжтой арга нь ван Дриестийн унтраах арга юм. Параметр C_s дээр үйлчлэх унтраах функцыг бичвэл:

Энд у+ нь зунгааралттай хананы нэгжинд ( ) хананаас холдох зай, А⁺ нь ойролцоогоор 25 гэж авдаг тогтмол утга юм. Харин C_s0-ийн утга нь дундаж хурдны утгын градиенттай турбулентлэг чанараас хамаарч ихэвчлэн 0.1 байна.

Ингээд хүчинтэй буюу нийлбэр зунгааралтын хэмжээг дараах байдлаар тодорхойлно.

Тогтворшилтын нэмэлт нөхцөл

Нягт болон температурын өөрчлөлт ихтэй урсгалын хувьд Смагоринскийн төхөлтөөр хийгдэж буй Том хуйлралтын симуялци нь тогторжилтгүй үр дүнг өгөх магадлалтай. Тиймээс Ричардсоны тоо болон Прандтлын тооны харьцаанаас хамааруулж турбулент зунгааралтыг тооцох нь тооцооллыг тогтвортой болгодог.

Үүнд Ri нь Ричардсоны тоо байх ба шингэний урсгалын зонхилох хүчийг тодорхойлдог, Pr нь зунгааралтын тархалтын утгуудын харьцааг илэрхийлэгч Прандтлын тоо юм. Ричардсоны тоо нь нягтын өөрчлөлтийг (хөвөх чанар) урсгалын хурдны өөрчлөлтөд харьцуулсан харьцаа ба дараах байдлаар илэрхийлэгдэнэ.

Үүнд g’ нь бууруулсан хүндийн хүчний хурдатгал ба үүнийг нягтуудын ялгаварыг жишиг нягтад харьцуулж хүндийн хүчний хурдатгалын жинхэнэ утгаар үржүүлж тодорхойлно. Дээрх илэрхийллээс харвал нягтын градиент нь хурдны градиентаас давж ирвэл урсгал потенциал энергээр хөдлөж түүнд зонхилох хүч нь хөвүүлэх хүч байна. Харин хурдны градиент илүү байвал урсгалд хөвүүлэх хүчний нөлөө асар бага байна гэж үзнэ.

Жишээ нь нягтын үелэл ихтэй хоёр хэмжээст хөндлөн урсгалын хувьд Ричардсоны тоо дараах байдлаар тодорхойлогдоно.

Прантдлын тоог нэгж гэж үзвэл Ричардсоны тоо үүнээс бага байх тохиолдолд Смагоринскийн дөхөлтөөр турбулент зунгааралтыг тооцоолж олно. Эсрэг тохиолдолд турбулент зунгааралтыг эс тооцно.

Проф Хосояамада-гийн боловсруулсан Смагоринскийн загварыг гүйцэтгэх дэд програмыг авч үзье.

Төгсөв.

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43

subroutine LEScal include 'paramg2.h' parameter(Pr=1. ) ! Prantl number parameter(Cs=0.2 ) parameter(delta=(dy*dx*dz)**(1./3.)) dimension uu(nx,ny), vv(nx,ny) uu=0. ; vv=0. do i=2 ,nx-1 do j=iz(i),ny-1 vv(i,j)=(v(i,j)+v(i+1,j)+v(i,j-1)+v(i+1,j-1))/4. uu(i,j)=(u(i,j)+u(i-1,j)+u(i,j+1)+u(i-1,j+1))/4. enddo enddo do i=2 ,nx-1 do j=iz(i),ny-1 s11=( u(i,j)- u(i-1,j ))/dx+( u(i,j)- u(i-1,j ))/dx ! Sij : rate of strain tensor s12=(uu(i,j)-uu(i ,j-1))/dy+(vv(i,j)-vv(i-1,j ))/dx ! s21=s12 ! symmetry s22=( v(i,j)- v(i ,j-1))/dy+( v(i,j)- v(i ,j-1))/dy ! sijsij=s11*s11+s12*s12+s21*s21+s22*s22 ! sum for i,j Ricd = ( -rowb*(ro(i,j+1)-ro(i,j))/dy * g / rowb ) / (2.*sijsij) if(Ricd.le.Pr) then tvis(i,j)=(Cs*delta)**2 1 *sqrt(2.*sijsij*(1.-Ricd/Pr))+viscm else tvis(i,j)=viscm ! moleculer endif enddo enddo return end