From 9e0fc22ba372c8e707168b96e627f457ffb4b572 Mon Sep 17 00:00:00 2001 From: minoue-xx <790354+minoue-xx@users.noreply.github.com> Date: Mon, 27 Apr 2020 17:11:45 +0900 Subject: [PATCH] Update doc part 1 and added solvers for velocity using sine transform. --- docs_part1/vanilaCavityFlow_JP.md | 4 ++-- docs_part1/vanilaCavityFlow_JP.mlx | Bin 2354981 -> 2354994 bytes functions/getIntermediateU_xyRK3_dst.m | 30 +++++++++++++++++++++++ functions/getIntermediateV_xyRK3_dst.m | 32 +++++++++++++++++++++++++ functions/updateVelocityField_RK3.m | 6 +++-- 5 files changed, 68 insertions(+), 4 deletions(-) create mode 100644 functions/getIntermediateU_xyRK3_dst.m create mode 100644 functions/getIntermediateV_xyRK3_dst.m diff --git a/docs_part1/vanilaCavityFlow_JP.md b/docs_part1/vanilaCavityFlow_JP.md index fd6be6e..8bab354 100644 --- a/docs_part1/vanilaCavityFlow_JP.md +++ b/docs_part1/vanilaCavityFlow_JP.md @@ -115,7 +115,7 @@ CFD に関わる方にとっては基本的なことかと思いますが、自 この手法自体は Harlow and Welch (1965) [2] や Chorin (1968) [3] などで提案されたと紹介されることが多く、MAC (Marker and Cell method) 法とも呼ばれます。圧力項を残したバージョンは SMAC (Simplified MAC) とか。この辺の歴史や呼び名は正直はっきり理解できていないのですが、詳しい方コメントください!時間発展の高精度化・安定化を図った Kim and Moin (1985) [4] の手法も実用的です。 -## 圧力の境界条件 +## 仮の速度場の境界条件 ちなみにこの手法は @@ -131,7 +131,7 @@ CFD に関わる方にとっては基本的なことかと思いますが、自 -個人的には Perot (1993) [5] や Chang (2002) [7] で紹介されている、離散化した Navier-Stokes 方程式の LU 分解をベースにした議論が気に入っているので、これをここでは採用します。**圧力の境界条件は不要!**と鮮やかな見解を提示しています。 +個人的には Perot (1993) [5] や Chang (2002) [7] で紹介されている、離散化した Navier-Stokes 方程式の LU 分解をベースにした議論が気に入っているので、これをここでは採用します。**仮の速度場の境界条件なんて考える必要もない!**と鮮やかな見解を提示しています。 # 離散化手法 diff --git a/docs_part1/vanilaCavityFlow_JP.mlx b/docs_part1/vanilaCavityFlow_JP.mlx index 8a25f3e36f0d115ebe082d3a0e934c6136768e53..d3f9a4291dffe55d2cf47a253d498a0e63bdff4a 100644 GIT binary patch delta 11393 zcmZA71xy}6*C=4z-QC@_Xwl;C?q1xX$U>3g4lnNR?oKH##ogV#IG69s&As{mNuFdj zvnO*hJIN$Fr)VP8hwf4TU6|C<|3B{q!Py`j4jJwoKkFAmNqBhIh;HrwYSsbyT5GtO^cVZA;a z@e@U$8GuRhMg|uB@_MQGn2BOD&)I8j{x1e2%hfRNk6^HXM=PczSY%wxeUX{m@sct> z!Iay>>2AwxowC+nJX7z`NwK_dx}`dx zf_m$w<;MA`-~Lsm^>Kxg)GZ90peKd>bsrR)a)_?|Jgn9`w6FCKjvFlIbo~RH*itBd z4(MR9Qk}XTSPWBjCq`{uBPW|>VHLnseLW+jo-6m>VDcE4Vv9e!A8X7B-EimMCl3@L zLb9<3V5^1FRW$5eh!1xyH57z0b)0)0ENi9)FjLn6?3_7WBjddAJnrOnw2)Dut{S|G!%{1liqeq=X$xgHa8`b>VjDt zaPmV! z1$!ibOha$&V5doH5O9o^oU=15&B;JyUWrSS+`$3uh2lCcbU8u0J)&bA6SMwIItBNo z0E{rz&+>qrduQDfBZ5->PK$?A+}WIQuFl^gkc=zDZ>WQl^7gJ7p3!)wPOsM*?9UUd zbMQ6)lUT;f7P_M}{G&eke(3=skj9~s(Lt0~I)VMjepcgdYwzE?O*3B+$}i2fBI-wK zKd;`>ugLg;PytTun~r|Ve8nby-OS5Yow^IHLYa;f24gy zbIS#vJN4(Il*y-E=7Y9j>RsD^neqvJTmF^4uFzr4yEZ3u=V_IvA@&6h0Jj1i7Ni*s z)%P$TBW4198w5(UB;tv$(~K-^rWM)jOIJ^@AdyD4OUhAlu^VQTx#z~2?b%KeG=G>w zWt#A*B22a8Zx9XHIsD%@c@q0~~A`Ve&UhCf+#xCMmJg3`Qr{imSe`p!=j` zlZM2fQVRFzMA`42|JkhqvI!=l?=4{u0-_wL@wYQ5t&~WfIK`Dotyg3r1qzyTQ*WyB zUh5jt5A0`FbleYgRe_|}&f)~BE!VBK@R#G?cKq)McTg42q7eC(od{{Khq2x>b+hg0R zLh^dZ82l?MWR#>`TM_O5Niz^HZ=oIK2ojlqyQ zR-iT90f=a#c3gmi(a!g0K88}VqueA)w5W&MY#NvdS^rX02`ZF!FAVpMM?srzuo+w` zR0x!ZenPM^cz#WA`ZUOmz<-w_v=AFo+_4y7YZvCNPhEzCh%GCx+sfo81F0;t0naw5 zZt-E>IRU6VpgU(Y)fdhYcCQj5gSGq^PnBEx#7zh$kxAebe|&>r!p(puXD^?*5;qVA zZ7z5+l0C`oGWtXg6k($Nblwm#2#)7;F`~WCgFa;Ub`mBOln=LhoxF4Z{Y-#)3)&s_ z_sQCv;}biX#5GSOqM_?dVLd;YlW=|m|&usYE z@gwWBa4vw1E4KHC3p)0!3iKminA4y2Vr&rHc>^Cc4rCQ%fxpYT#7_s{`WwJXd`Gt(Y9j+C?~Xzked6Y_}< zAzH9NyJos*>8oQy%4Q4mz1_hv9Xot%v@^;h-HL zcRV7&wnbaeRiJ-@z5mb~a^ZPjt{3?)3qz?*O*gN!W}AtwoYx{0+$bA&9P4orTmNa} z-Xog~j}J#JH^v)Rt#gHhRR~Sr6}IfrDuIYZ zbR;}rIn7}@;%1;lxgykkf7pJ>ntCb(k?SR|f8Dl8P5W9;+Z#IYctY@pMK1B8E&|Vh zxsZoOstfuQ?}f(>uokQ@ptR%OStETJUs<0&m)kFBg`LR9lpyM09LOMdFa*ox!?bOKyH4^p0Ng8DE~ z$cVkd@hxSK8yKB88d3VOZjvoQ(Ais%>H}$D3)_88eiGyb~ynxQvQ_G)qj z&rZl3n7`uZ`kLE-UEAjmuj|kUMJI_YnAL6cs^(7bMt~3wwpG?Z^;Wpe$AB*``nrqab;ZM{r=w;Pg$(%3JTHxpm6}c!qel;DyrS zmt5=jyCw{?>yW%KCtM>Y4QGs6V9A1@sex{_&I`>u{_}OUfk{3>nJM9wa%gqZ^>tKF zLghRN=X3USmFF6(KkZ1&NLADQmO{hVI_2)U7sv(Op_ah~mCuA)TT;{_HG&8c^tJ_l zIs=cPU~f}vm>A1d+DS5ibA=$?(qhOLir1%xGo-GqsFj0%CXWTu{@D8R&0&A!ua2)(z?fpVd!3|3{T1&;y@$j+yNJW1jSVJ%_r6!y`#un4CN-LxHfJhS3I*GG||&WFi3nD>l(3q)b|q@<$Rl38N?UCY7y9t z+tR{k48>o%*as?r7|T4Nxh|*Q{v6q`6RMgp)mLm8UDe96gKU=zI^mX*x_aB*V?))I zf*s5iY$O$U(t?4&a|w79;b;H$8Y;zTEN6Dl>FTMhF3u|E z90hJNreGQdKx(S6x!NO1H*4%OhwPe;abw%}HqpoHmxSwc^?m*Qd@9!#j)&i0!b{80 zm1n$vTZXkPKfRTnakQ-ltC>S=#Fn+4VJ^*s6QR~IkucLHQ`tz3)GcX7paxbqhaxgcq6kB zyN`rG!JHUxq+ia}ae-8ga-i;Md8KS@szHOo%~CDeE{lP~n#>?`@)jyHA6 zzY27E0I1+Gk}lFkuzem?>oQ;6r}tRqJvdb^Xq=PFrha~OPH8|2as**UT>Smlu9gda|sxw0@Vlz!R4{+LpISIrqTnRpY| zbN3h;&B}s%eZf1J-}A4LQSow<<*$o`+Bjaj1#q8-A>Yj}?it@sAir%@u{(b#dZc=b z$+v2Z)-TEF#O?I8^iwgNRjyeu!F&uc1zjBjsGYC16@Q;96r(1?%<)RDmdkT`5&+`E1LdC{|xGPQ~Mg80**@;cM+c*#1Xz%%-hd^Y9= z9{nShV6g7j)S`7yn`ET5e4J_S$eYsrQp4vlZY>Kt+G$ak$p%p=NMF)Ej1f8c$5db( z*iu}QR#c@S_@Gt!G4R-8V=@kkbV~Q506ON8f|n3(uP>R@n(~L!s~PbWIB#m{3fFCd zQB;WDqM3_CaE^p)fkwniV^2`K(wSpOh7a~I)x;5bqoK);baNsSR4UFo^*(tv6K8|u zpWRR|7PW%4?|DXP!-hSXfXKDzEnBcKlx^bKlTM!1@g)`B|9PSi})8M>kWOVIQ7uST(-qyQGLF8o zm)9AYLrD6JktQ)6&d6NYZM8S(1JhG6g-8WaQK(#|EjF0n%!&k}6x7QQ#Db}fy64IiJ3m^wuaBa9gx z59F99f$sn6?Nl#e#3=TZdgQ0hM8)SSF$r z1UC~q%(fy>(z?(0ENN}dn2emrJH!YAJTpMT(-G_P=C5rUq;UdUxp0Xtj_HLscgd@w zRjv!OB4?btV5KfbpW0YrfsC)V86SBp+jJ?bZLmTj!BPFNXL|r;qYcR z{c@L1U-8Trdo-SW_AQH}Btr|h%M%?vBh#L=A?9&hJ+Ts0P+9U6FnYHQO3IJ9;O+XL z=YeFMnmXF0uZ-jnuERdW)43#?z_LkS*T_AQ;N+ft67Cl4*D*sMbPJFx4sDQ&jcz2R znd!&Q#+|Gc3^cYu_0q{~g2>_V#l)owzR6?hIT`Dg&BIr__fGOq&u(bQdLDMfdtltj zl=6--7wsB~xVS?GnB=#XhKyPyO@4m$`IqEndEXC(TohN5U&;L80v1#`R8r3xS{EiI z>W+l|9z-r?CSxchPKrN*YSf=3TR=0oXW)0h4}?Z7YW@;tuHzcL8VupFlp;Oxb1fCe zKG4QpkAyXqx2#`|%WMs04xacdrULY|uqG`(^p%Lpb5YxoQJ!ojoZ zSFI3npRu{tU!e+fx2e?k`%Gf*cQDvR6-0<6dU{`V2l8nbLR_1eCf|!;{me+)Sp4Be z#J$rW)&NnZLj#Obf6-D~2`GV|rJ=6+cE`w)gqE>)KiEnUpXQTt`mV5p$A=xyg-z9P zMUEXSc&)ku*r*1dC*th(a^ch#ibY7sNUXMB{5fr}ZO4No0%D%XuAxUAh)aG?Gt!a$ zwzA;(Ybo7e$;~WpKO*EKVnaPx#w<}7HLEA8B3jjo((0*2Iz*8acHYHqnnY=NO@(*Du>n;s35#wuB3y+(%}Ch?jAYQ=S$j93aNfk|cwrlQ`sz+!y2w2eq*a3b zQT59$r)MZNl#KWH#cK^t!;5x~5x;as8ED1@w#4Nm1)+Y8Cm;RF4~|u!`yEt&Yp2rn zTkIp1&Gv{@_1nv@mY?=C`dZcK-*t1ae_YTONK}u%CX0x;3Z}h~+pnnL8WR0w9>0&_-cTV~^;%VnwV(KN2Z?XD6IrL**Y>wyN)|R; zkbhF_#*V~qj&&tzjyQ0t!HCmJq^PaG@csM^*awVZTHc8)5mw;PYh0uHOq&yz4qaR=52F|oFKX@12O;N(wPzVgd>C5omjmv)i2T|GUZ&WfH zLQ!7*tp`GgmY!g37DyMx*XcTsF=o1QJ>N0Se}(=}Gn(;hV^#w>*=Hz-$+|(I+Caz) zbDbu|(ODLixOUCVO^x$R5FBXT1y7d+8T%I*zDpOJd*p;^KPvH0Yh(Jv!zijMl`j<% z<&N-l+KI+&c73+s*iwK4MHx0D@m5!a) zk^Lvo^Uke3{Y63IPSh+6W1CZbCt%S@Rfio;peJL;m_x5?HY6lU9DGpC{^T71yWE0 z{A3uEg1Lf|WaU<4^Ap+DT8f{=wU>9MmZ+Ptevgm>I9CDVZ;oC0zKZ?Kb)a_QWFu?d zexMOTl1vBQ1@kxTXzz}{@tc)kTO{J67c*jw-jwV%dDqSFq#erht6qA9(!yEa=V-`; zvX;>~+4n0+JV)cfHWl^p!t{_H&xjcW?s1pUq(S=je2GH!AQY~}%8EgPny@4itj=A? zBGt6VhcZ*R^yB`)F?7Q^Lm-};l$@&6l<_&dvx%^^0UPIs+Rso;5}aUIC6&;IS6Q72 zgT?E(v}2O7KJcThk+m0Y*AEqsVjpjK&TrwV3fLApk(h}p_ewI6Lv}cju2uLeU1=MP zX@s7L5JL*_mf`herV{BH((ag+9JA+z|JEmjIxL9A{hzCg$cZD#)&Q*=1SMa57g z!A@&o;u_HqT+Ug z4$a1GLKu#xW!pDIuxe!hdtJ0(V(<%98LoWG3)h4JJM*9LBoKC^DTb!Lmgfkm8FF7` zTfd8Mjr0IJi%0~2e|(}kESBmdl};2;CRVTm?& zQd#_yrZy_6_Xn>qYlG3APU6VpfB$M`lMU>NdRbQIOlJ#2UUt*WK_55t&-L0CxRnnS zIwN%Q55JPzt^dAqj>)wklf5^*zQnCTyP6ljizcw#_u>CXP&Y&FzApFKvFHTry`8jO zaGsn+$gJeZ-b>n4PUIe{Lhy0rhtF5;FNW<4H0x+pDYcjjp`)$%_nj$p!F3eB{j&I~ zQ&tEi0R*xI3D)uva;$Ik!92&`$o^)7^*y9}5f8^JyGA%ePsl(^ZXZniWT6nIN3 z67`!&omWSx{wXhgL{qqmcMNIreKqgz<$_@>w5z9zO{C>FWPdWxN(JVfOHVjG)Mp_? zY6C9FcIa3<8BI_tJa*eUUgxIL83G`<%Bwwy0Ks`fH9f~owc9q>_7gTo&*9bmyK-pX zPA%k&MYEE8Vef*VALVT25eOV)J#lvZh*gj5>7vD9C`_z(N9%rT7q!=8_xJ(_5P^x6 zdZJ*f^b*%>Jtrum3geyW8Ma52&y^nX6=evJ0q<2%jmGz79%w#FCe7>@mthdG5K87E z0QpK9Y${oug(h)zS%mf~Z(X#RL`njmoZa1&|Fb>RV z#yM?8?NiD+JHCn4ITsX9hKCVYR--kb0k!3Dmwn?f9>XN#dUsk6X!Iv0YINu~NaWba z3G`k`EW={cwv6lK)Q$BeHcADMT5*f+PhHW486i^w%cSu;jyjZ1qhetf+AzfFwg@>br zWvhB5yEASJ{b+(gVTN$?LE1JwK>i>jSB33yaSz4DzbyD>n`z2u|9-X*&L9Lev;y_U z)Ue4b#Kq>@1~~Nm^JN1KGuF&Xl;*H^9o*94UA3-r1w{Hgo6uroG=2)6;&PuKCGr|p-;%FJBDGy9`>?Xna#6$rK)qD ztk7LaN~LA>rK#QqVm^LU{dgs2Xk+P(z~+#^G5n1{r>%~{_DpLrz^T2Yt_hgjBG}IZNZ&JU$D2 zEXEc=fBSfJe)yjZA};30LP-2C5Fbi{#Wp3Rt6$R#oX?O=4OXz_ zAg(M5{X%8e7DB>ncyn#Vy4+)xz$98zNBCi13a#8sM-GFgN(lr+qle_N|Frch>F07^ zb9SK#qkmk_q`65HKiIvhLtG&vz*dCQ)4k-E_n&hf#*1m#bE4?{fa^m+zsy}Od+e?k z_(3Fu5NgI$OQgt7lRCR>X_Tyj8X2<<^;kaw01wF)+Nl5PG$^)adUmtz>LB!C-Y; z^vy}RuZ_3i*xIo#ij~Oj{E890=Hvyvc=zubZ3t|pnJSv{u~9jJ%QCpZ304!ec1dVb zvvF33jpLip)W-Tt`r71=e?m&L@6V6|OlqK#kFVIGq~Rsn`JF&;rp>?E>tBPr(WrYR zmkEdt=>SNP&JfQmGo@uZJ~VOJ3SjCvix2e9;?OafVNx86rDG~+eD}e?6+ZxflJxoq zjrsL$9v_K=P5*iG`**g|r_>u(fvxPwu~{fIl<_6s3wV` z=&O-gWd7~^8UKRLqgAbL!qS^DXwoebP&)^iUJ3*){Yo7Li8;OZ9rJnI*|{f1L*UzG zd+1$fG*C`)ERvu*MC@Y!{=6IHWzu0u-`bYjIo8-5fTUaaq4F1AUqVnZbGeh3?u!&x za{_;5B>fzR<~dJry)Gc*U8%p&E09Q`QPrcJUS4bDK}{FnrJY#B1XaeCs!LK8MHrV8 z5(r4SkMEd3c(wNt8-l53`l82LWz(>`Ag}x>NUE5#vYD^&{Acl8bL)y00YasdJibQ4 z4`**VkaZh3&x$MszaVK}PaL8k@dwUs%)NgzPW%KFW})Sy!V^uu*UYKEKqr3%tNMg(Jv zcDsUU=_#VV;?4ckkTopCB%>DIzbj3+68Q0omrXhx^WBhnn&_c^ccb>iWQgCl`U_Iz z;P2f&DXehfDe40C-Ihhcek{g}k$S)qCOr23hRHg0gX52N$qne`22`Wwa8+nN*$NZF z?$J5Wkt8Bgawls5S8%PF-5LT*q)>}#L~8Y%46Rx(PajN5HZ&5H1g5FJa6k${V38Hs zN58uX1uk68U^spGzrjswT^pG z0ppwg9r?up>`p{;5rRa&b^xhruL3NxgtJ#}kUm!1$cO9KP-f)Zy_?`j)i=RgDxV4w zu#k5d$F^TOuXEqB-GsWx8{pa@%Z2X-&~sKaqYx{qg(bPJ#oJkoX4HTygWh?6zLtNB z%8MXJE;d-P2Z4QVdl5?%xGsR&#gzdRpw^UCe2vogHov#1j zLX6F3{^-6Gx}TjiNk&`=y%NX|kc41FSAnQ-LSmuI{l2mpYl@z2sjF-yiBhOVN=h^G zMTT0PM5Wh`w#C2u15XtDD>9jpdmPVk>w-+B9A&h12!yVvLdN(Rw=JuGxr8xRW>z6z zs@gk(l3ErAQXAVYa)R(eDryYtmkHW9d~X_LV@y-+S@o7m)|SL;S{}$)&w)fnIu$&i z*N-o)!M8xg?ya#FX2zL6A^lFKMI4lvL`RNmC>@JTxW!3F=j8<_W#Xgf?1=tuyjPVt zPR9l0yFpP}1X(G>eaAG)Yj6zpu9(M$v&IFRX30@PL@z7f2v00}0bBav18EKSLU4W3 zYZyuz-AwE?dgQ6zgV;dFfHb+UCEjL^Ti4RctNHq`8hfiPF~&#~M-82_fK&H3)`Nf5ml%Sd4R;3MD%4WXA6ntzKoTW5Z z942Qun{iJ|{`b=cnw=ThH-M^mxQUJ<+IpYRBl5M!i>??PA+nX#ziNdXVMfrqpc_^{axV&;#$B^xtVT~$5 z$Jucgx93Y5mOe_wKwtWVw-vHLGq^6SN*dYHWi{ z-QODikv{nT)1qkX1jmzewCqFgi&GBWds%+lQ=#O7IO2C6oOcqmL#s?AB`Uh z!gr8YfEB#r!(IXEXTcBFz9BDW)PF9s<6?PXaG@dLnzk1TJOxC-@)23+Hh~}EeafK( zX7+9uVO10RL(v%4HZJ8tO5AQ{t{2O4%t|?roE3wZdR`pY*oDs5I=P;bAm}V|g<7{Y$%&%Y^*G~X5$KgE7GF|sZ>w>4_8KV%;~J>jir4n$CP^wh=~ zIjWAUb@r>?qKCeybj8}-b_6!TR6=Y)q^%rbd5oc^y;z!PS!x4N*%$quhzOd#_Hiwq zH!wAevCY`N-oIB8dbt?KTDGZCqvDpwAG`EjR8mZP%(m;)(E^bqO7zh+;ODa8ZG?Z0 zTZN;X6^W!0T7%L$#8ucX?S{2;h8CEpQ|dDbSHnF>3$oKn9hnem2!TP81qN-JDCxHF z?EMfDuDv8WC7J`!X*P8Q!GcoCe6RWawivDa!6Lf~6>bosY8!Q2%z^V>gWPi3D$Gc| zB7-m&ZAHDxlg4{Mkr;%K8U0^*6?VwZt?kj3pQd1+g*_WJ>et`-IR_4BSaGK@bOQqN zG3;7jWx`24H5?8T?y3{ znjPNja)$R^MDZg5j~G|V%Nwp?Ncy-Qs0@KS#4FDtiR>53fg(6Yl52-eo)^**60wyB z>+FYJkx~q?uMfomjdJ>|Yy*a}uLHIAh`ZpSwK>q#fLT>O2K;$96&_?2t$L^}cxdqO zDAu@NUefw5XM1h~J-|88`5h;GO=#|rZCb)=4qliV1ZflPM<;aDJbISW7e;~P_2@7KKxNe$y(c8z zrhIBL=JCyN-o7#uabM1bsJ*l4NCc;JSWNf&riy=GRl4I)q^aS6%(isXXm!huth4HB z=G?{{XueH^v_A_&wLHJ+H~#4|IPT3cCK&x!@6SakY~$)~>3Rd{kY8?I{=mlSe{W=E z|1tniBIWn6jd1%mGsX0LzKTIkKNxYx^d3x!!N%>Dd~%yjNbTz(VAIdgoC9<5S;Y{U!*|IC z)DG|0py;fooWLgJAcnE`)ps6&mkh?z9|G1%VO35|vA?fQ1OK-Ubtc53nK^a}uA?9=d>mxe6Stsd4A^TyMMV>$w7XC z`UL$6<`e8ExKHq(5I!M(Li&XK3FQ;&C$vxKpD;dQe!}{M{R!t2?kBuY_@4+q5q={2 zMEr^56X_?iPvk9aaumZy|H+~yn2QpQ{=dCg1Ijsz|Mn~gD6c4?jib#4aX@;9l(azn znBo8%ngl?7BUE>*OsAQ^B4Q92{0%KE%R)&la;eTve5GdySwUA6mr3g1_+P>Dd&B#` zxZ1R^1_fN?bPyA;wMS*tdC x3YpcchRQaT*lize0TkuGY&^c!u_R^*KDia!Qi zOdc%sia;g}y-mHhseo&q`EnTUg1i8~i^X8ANtLNzhVvL}?BApS1%Kol(mC&#GCF%NZrtf5*zs zC<(!}K!azLWWeXU%BgL=8d*oAU}T77U}&G4<^Q)PIY=;j3v(M2RvUX0D+_iOD;rB? zIdBLtu>ZHIlw%DbsMhJgz!=%V{#SGW7v=wqc5<-#U(rO$I~FSTWfovy$&UX!JQ5Ue zPDx7opA6?T!5I4}!N5#D=lP!*nu5B|DJcp6yGg)52R~CI{J;5sK7KY8(Cay+G(_$Z pDE@+y5QK9cQ*HQRZm)OQgMdw`ki_dF4WZvYT$Fpl$aQ!P1u3(Melhr6|VbLta<2N-U-_>4kk zU&;i|bTXGGGdq(nMX0zDqqMJ)lFTzP^M6x$JHw}3sPfuidgz^INu0PJZ_W=};?KQF z=p{rBYi{jcrO&jSt$O!os&O6m})|HZH0L(Klp$omNUH_FjpS9{Zx zVxdv*i3h z7BgxiuhP!iKg|4~W^g9^E?a(LTE@sT^WQw|Rk}}MFv*CEwoQrc{ zcqPxw-Oqwcj*O53P=mIgEm)bD+6SDZM3(XgiHRRe! zQZLUVNgT^hasx{jVu2kbTCjZj5`smbFR3b(wM#lVMhggBDD1t(zIGCYXRXA0Q3(c4 zg1_RjxHjuVW7%roz=*Q|66KNyX=Zt|c~_sH|M}NoH9KGVBY}t4Z+WK2_{;Iglbe?+ zdN4wSjgovyz+L7OajKaj_wlFN!87rAoj0WrSx%%7dD=oSKgz?rqkzY&`GFHqh;-}h z#+T#$)bpT?(l79qd0F>c<<#&D1@~JPLddGkvPyRKqDkQ0D#f&6rb5F32(pkI+?|iPu={w$dHB)z5c9-nd z)sP=&+TQ+%62Ia<-rQMRW83h7@onDD4E!7ExH%AP-KnM_4ski(z4oRNMf?h^=#kE1 zR)dVk+QM#zr>wKuq-`t~k8R-Im3MBlWRe1=0}tId6y2uEcIreF>#K+<$dDYj;knRt zy-FzXyIHm4LM_ZB^2_VKlQmS!5+5z46McirpxQ=QVf{uAA%J?1y&xX|%v5}TG~kjw zC-HVwM26TH3$WKso8O=M3bSDtwUcWoe@_NVHg``VARsAYLcxYLE+LA?Q@G8 ze$H{3W!<6itDMq&A`XWG7Psegv(@!@1;W$8^3m+WLKeSwl}82^CE;I4qD2NXJSp)n zEM?;y1z@WX>lDTt)LZ?Nlz!kL3xTf21RJ!H*liWN>kDs@qYhv|ERw@}WHbBN($2bU zMdw+pztO@Fu-2l7`#drok3)B*vY{+hZx7mkb^wSxudZD;*%2nK{%&gdBFZ7>Pcqlu zoDaEJS3w;R#%>Rz4j?BI@$0te;E@+*w;dB|b;a|B06b;iiSsBxu=(WAb8rYN3e0?I zw~wTguFR9JoTyJA_&da^Bk2!&In}F?HSban)s>`#-`In?G zPP0&ncD+hG%0uyw38Q8KzkFO~(e${QNAJ)@^`I7tLvu(b2*|v7hqkRzhriWUMe`z` zc{IAjNl>7z5G~I_ya)@pSjlj$rsy;LkY2D%BF`+O_7w~RmZEtLR8zbSJZ3CkCdHIm z_=Ei|?(vtp?U{~L?WW>E;D-~5goXCmN1m-Z-hZxMDTSKinv0}f9=*q}aKufa0CJl; zfIGf!N$-hV3n`h4kG4oMI&>_t- z@J7<_RB*OA+;p7ZO0s`jE#u`vaMZ(9i#v?#lE`Rd!!w}9wP=&%WZMUh2>Qc&O8o=Z z@Lgre8>{6HOGZM(0cK_r%|tnqw<%SJx$TeGEn8Afhazyi7WS=Mw5V!b>uPyH1s=-^ zdD|5bF6rEG_p33vt7kZ)P4is1Zv)pt^!OFG+`4*%FB5ATixvv{1T3JEc^Q&K>+Fq!9R%yKAWA(AwH3Vj6G`^kC)x?yXf#Z)tp|Z=K@**E~c-?{BsjSvt zQY1GcIs$|8PO7Kq1>fa8!1zoVA<}DN%43a}IRZDA#q??-n z`$u>)4g&0^Xu;4JBhtr(+Y<{&jW6ze3Is_vBJi&&%-!9o)u0UT=$~oi@J9xZEhJ)^%VlH(P0-TDIY^#Jk@E>^EuGB1vy2q4_D_Ug4CQt1O^^2QI1Hfeu)(0k`TTvQqz=ERm9 z65B0;SH)jkgm)rUCsRq@+mz^C>e;c_bF?<$I@RF zqZeorK@Fy~%y!op`HBU6SX{tGo@p@1Py<2gL>OnMqWEYZZ`$r*KOH6PZhaQb)vsBe z;)por7!&yudhk<=O}!65f>-!iTUU0|cc<@qhhs5<)aX*#jQ;8`;-h;ca$Z<%S=-^> zZ*n|K_Mal|(N1$C9n$iCy(`22#E;o=e$M10lP1Y+?ju8B`12VYj}j%MkD78uTma1o z`f*;PS0(CSBz8~hi8U>_e+6sFg=jf5xX_D$|80e+at`zD>jMYcnzfzr82!Uf3Ya{P zuIL84bnCxM5su`Bd1@I0qHT&AK&A#%g$7U5K4y9^qXV5QH%R8Qj_IkvMZ$E`Uq}lpg?p3cyn7ZC(AWKlXU;qP_|Yz z{@xeJGL}aIG>D8OwjgO=N8k({SXo^xL}{yaN2#NpEvPIANh=AEcnK2X$EA1O)KTI| zl%*O)4>DgkGiA=lWx1Ld2s+-`ptx;sK+egqOFNrH(l0&Nz2e0CCPK2V`DD!oqdQs{ zU-R#7U2$I$97o5zn_1BOLye2<7w3(geEkaKx!ZaeijmjT2lXb?Fi>1$eBu(DxiMdW zG}gNVOk%psnlDRdg* z-m`1WPo5YT^=FzbcY{6p+ZHW!Q1?O7{<~6sn~M4&twcvdmFs4ajmc!8?B|rq`O5a9 z53D0ck{I20*jAf+1>i|#FqKDTYj&hFcO*zi+9TG;_K^xS&>C_ab5$dEyTafYyT9~q ztSN!Odk%yDS=l%&)S2RsMjlLa-8hAIoa^5dTIN!8cAP*{^2oH$bg#QMvv%5~ZsAe0 zwi4;s9dGng{^BtU1L`| z$;?bRZre!{xa5%e2qFdnMr>kj=qq;8%eXYh!4M21e5aL=Napjoahnx=A8hDk4Bd(- zxa+$N_+7&{8(?BVG+!NK()8(V$EUZhAW-pa8o{VR!s?;i*aS0fl(Q&Ej(?2$$D5Ef z#xRkfLFxMe9%B)!NK_ydiom?d^d=K>Q7HkbOx)4Ald#?egY#=o2;@@hH3VJQ_Vwi z^IGJ$k0OI*bxh{L<+!fRbYvls9HuNG@s5kO0&IpWYM3J&Dx-z_K@lY4P)~hVO2K^6 zOSm8E?hrGTj7FQWDhbLuC@V3bU3DD4k}L=|(++of!l5te0j1=|jLOw+AP(h4=+&)Z zBHpv8JeqofJ4B533(I@%{AZW$vkz)1ijEn;_!I4b?Q>TMhIorrDjYLxE*?`H@q)$* zjB~(t(1suro|~)l8263I)^f0dy;NcWB`ze7hyzF81i`tD4FAIF$Y7r&Z#RY2utzbk zmpKvB$KpfcNs6>*8)38#SxHx=0&s^fhfS|pCvL$lPtpQ@a~I%G)*^brBmyP_`mw9HqtP$CMyWHRgFqhyb7`o zzejP0RLbJhisI!KEs}nE@e(~py!Ve}X|e9%rWMm)IU{h0i!u(3!N2;^t@#~Xi(Wj) zk6ui6aHO_c1-GNcf-#Gg-@%J%)mhlt9lw~nx3U=6!sQC;m1zjlCkGiQ@9=g2ckco% z))2b{go+g@8dgW&_B3)#Olw^@zU5oFYX(+fPb|pIeLgipSj38B%!4sy5{xo@=SuR* zkb2PYY~}7Rf(R=8zBOpC5pA<^C7+XiFdE+R)**Du7iliHc&NOp>@_vveq7vRv|JbA zLbRJ1@1;&a-!ac!zZPIZ3|yxJ?n6JEt*h>2rj8?Vq=@E4HMXu*iy$`Uc_w-B29@&O zSpST;B1S~RQ@kp|c{#GQ-$SR8!AqqoyEtKAdS*B6YUO*giY!$QV`l1NWq(k*%oii1 z6AJC~M5xT@xdG=bC61#M4urV?^v{?~Gm^@x9+xAi?E~&>_BBUc580oAs{3vI;R#oL zi&c}#{mAUp2%lG~$@rA?U55;GW9mjG1dGp^9pk}|NE3hNzD1xU)7M}Zp{~h2{+=!V zR)KaXnX_*-$yTb$)BwtwgBFm_;G(pmEHW^94qnEn!yCekyj#1qopwMU9;u8dp=H(# zMX7Y$%dLyKL+N{BoS@_4N8CJb~}mC@6ItVq&1~ zG_h$~vX^K#Nnn+=n&t@=G$-{fq=Bmp9oH6A6lmzdtaVl^8d^yYJapr*jwv^ld|R|b zx5Yf?NoMu@cGzY*ERFG}L@_`3+yjPh(dU;Ty5_z!Xi3$W#Wz76_`2!=Mracr9j8Ay za4HtZ)F*kEQtV<|Is0toYHW*HJ7AgcEP59;6bl}pmF3_S=>v9H6M^Jut#sTKVbTb5 zZ7og-{V76I98F^~^9>yc+4NBLs{k@d_m)Pc6Jg#27?AFzFAr`ywCiX@szJ{ThIRd# z^ysQ`PfJPqS^nt;yf|tgw}{U{(Bwy<6I8729~uysH`&!W2e%@}&cKMcK;y@eEX=mB!Qg|>XG!Q(zT^f;M-@DE=Sq}Bi)D}{gSOIUy zun1F?Rm5~+hH*?S5g?o4HH~w~^F~xoZ|egl*lao=<~X7UjP$3Gd0Gy+;iKW7#!?fR zc#MZ@3j1J_ZU&mtFPrzd?c)`e-D?*f`m^bnIvvy@(x#{nBS*_2EHHj*LBXZB)# zXuDjO-wuxmZrUH1{D>((g3k!MRvU&HKgdN-<(lmcIN-F>mGxI^{R6p#QsTw9$)e?D zvkBjV8baO$yzN%&?vQmMaFxQ+Vp8Vrv{viggqKwoq@|fiU{{z;5Fc%mm@rw<4rDo! z1RXK{#A8^$O`(9m6;`QepBn3hCCVkKVb?Y44)RRpJ`VE4QtFkzO>Ab?_bHop4#)5i zrEVQ){DrHSsP{%to~5HWx4snl?2bQZu%DZ!60z0>q%_y#B2&f1X+rL-0P7#mNP0C` zMQ2QP(Ei>Oc438Rso!c&pNd>RfQC~>!Fhn17O1!jX)x+I@x+o$>Db)&+h>KW8(#~#YR zGLX5oKv}ksdpgEXyhH8?#}&r#!u>fVr{+Y`rD&zLWd^FZnb|PtG&V19sXjcjCY6$x zw16N^IOt$B2)na8O^dO?P?2NqJS+>`n*j%0Of5ZVr5>&{kbkCevL$dowEh7zmkuQG zDP!QXrW=M&xpivi-PSHxHY(gIiwA*JFw>ev;rWv@-;=tC{&49KV^vXv> z$i8mFiwz(YIQHo&yA&FmxMN^83Zp@5qK4A`yF*~A_MSYq!KhU6PL*!BYU`ioSOkO=gI7y#VUjJm80} zD9DlC4XsTjF6hthMy!=-Idt0ukbom2SmS}9Yyu?G?uve79GqI;9fW^8AER^VX!A(U45$Bnvv6iLknI; zWlN5>Yt8wjW%62cO47Z+S95v=1~0Ntcu6K#8g%k6yx}D&E@9LCgwgpOymfzNEIX{p%16TL%~KpM6WRv0?The<}!UvD>T8JA2HEkRR9RqtNZHi(67kjoJXTc=fDwK}G1^akcJeX>KzJoWV1NFq-Z?`_ zD=h<#%A1IWe5YU15Z_=j0SLdd?8KCx! zV9??W<2TZ@uq|9go?$ACfv!Vc`OuK8=D!mBqD6G~zEve$?&7>9^iUR!koG#(;YAgn zekrhljw2{e-Vm}`AhM`$xNagiUwNnIcC!Ec7I$oT&WRSe`vkgk)&h3fT%^&Tvx#7t zzqx~-_7rb1!!0CX^(V}M7%GCP%G3a)4e2f^Hr3?7gtq~K-?gVQgzKl5f|f{MjdV4h}Lb)PoE>>_vqy}m@uL8jCNo` z*5|5Vpxjnp6s!$QRD(caQ&rQ}!_8S$EP$g;8DRvCd~2A7M+|fP-H$k$E1Ihx&=^Iz zUIB3Ey5@`9z^tl>me_Cqc%oKawf_CEu6j%3nhDB#gs?ygg9y<5sVU2|b<;vq!+3^N z>0bwbsy$_0+xojrJmk@_-*r=Q6MtOq{kvafdXb4mM2efH^XpzYX=c`E*kv*cHM>!s zHiC7N-2OTxdQ>bvr<=l$1s0dkjL-meaZeO;u6nTtL*PbHPh+^02uikSWP$~3rDcTN z9}`!RwCT*ewY&9+aWrLlp;<_T7{oM- z{ZjSAkDl6l7)@1FA4iA{>uq5d&by9po$FQzAM~hIVNk8W%QehGAx_Tls@31ty7wS& zL+#kBG$8VGZLyO&p@|8`4P~>cxJmDHCEe3Ey41Dgwm(v&y8f}{zfOOPxwiUr@x(c2 zx3aE+sNqz+<^8>!^(vVC}FnGW#dcVioQWmXk&InH3eJgJ?B#e=!kaJHG+#oK54 zQH}QwZ72;bW+=ZECrN2|Ds>K;cEjO81q0?K2G)j+JBkX=!Fi^cAqm-Z8RNJC7RX?H=Q~4^=xF6-%lKEuT4k{OtXVc#0K2Zw zm&wxLl|XR7@R`MM;$eDQr7O7yL*CI*b@a5Q5KDNaCJ!Z-b8?vj^(sZ&WWC0I%R_% zof|l6=7?JMiE^#0j21n-8QmIaFd2fJ71)6&S69I}J91yH7o60B+};zZ0$qO=%F5Jb zC9f5MGPe6S?0RXr+22d+a$q?NbLv1gYOA2j#mQWCwz~`0fM1-a4|>` z9))MPgrs1uM~!o@5u#&()>^h+0g2Y`51+^U#*1Ix@6Ja~Q11POUtZb!J&e0n8IVcw z(}yZZkT2gT@+MH(xoll=kQghDTQDQLTb8Na>L*{Jqfv~L&W5NYo`PB4XTw{SCiH>Z zsA(u|6Az^Y1OEBmw(QO^ z=4L+xjfRi%A2>ZI0fn5EE*=UlNe&G;!oLF;*0|LV_@dgiJtYDfOg3A&(!j*3zXmk& z>dVZ!i0~vte@$r-M8;6c(j@73q4z9|wv-8*JTm}!_KxA$0T|Q<1G~o-a4~wpt=BaB zboaVC_s8Ggj{{HN2txYyA)=JM1^Ok5q6pxG1ZcpPLLA;s>AlMob3~#>P(Pp0{|LeV zsz#YGr*_!AfmlDBt@KL&sJA+{=9$z)(lO(4 zKyr8H!W|}Fkkb^QCk9+`H#JwhH}uJUtO$c8O7`rr&FL!0MDr{@iNK;s$$RccHWche z#g`?uRO9uCc0pDXc{1qL!74xD=c;V$zfMp6tL+O2h_G{q&~sugy~>S>RCnpp7&Nr^ z>cs6532&6l!)Xd$euo8Iq{o!+WNA^4Xn7QoEDy#?2pY6rnJ196&$#81sUe%QyTey+ z!tQFf)`180>4cAA>Od6!zuPV?Mp_J1Uq70Pq|Q9K&`R-o{0`8{v{vNOIN7CqYvL_L z`tQMAUVePu!B{Hdf>qQ_0+rah`SMeQh~iF7~)0_^!Lzt;e*^B z7Y|~Pu?ayNE_FRXQpvm=W$zi1^J0K=K_s@_8sI3f<`9!<4*>t3yy1+}^i&xsaO~?XW{&0xqXs+*Cl_p~Ma9<8 z1IAeHguCi@152nsw>kr=-<|C*>k5A{Zdy7bhn%_XbJ_`;Bg6itbatUj4^VB%Eo+F| zvYbly#9>OPL6S^1Z=obm;$rMEwi#}AP~N{Oju7Och?gIkPbnYo79tN!ji3p{EPcjn);~ln2HXMCCyz)bWgDyXu`{{4GtC zX!ZqTtSgiEe#qh!S|nqEpHi43n4t9tEw2=0MFg?zA|(lFWTM0|b4^kw5cx47873QR z&8xN52z4euQ*%Qo-t{K4(q*5U z;cv0i(0F>nO8Wcg8abeS81CgJFVJuRgIpoW&4MlD6Rf{6&#iFxcP*a4LpCD-EYqjU zLnGAGF9$~F`~bGUe}i#^_{Dm<88Y|BOZ-C?FuUWg+6>jkz=S@$Io-A`%3jvhBcT8A z@j&zZ&nZH}pc8{ zXKZvozJ9S-?oHoW$Y_@QEWJrq?VNyK8&+E$pWHEZL&41;!Il%fp_7-d9Z6xR3Gu$7 ziGrk?+fuO>R@{7|+7PgR4Zv1YDg#^1Oz$23eYRhwF&n-bGudN^91yT=pRT{ahm8Bk zv(oTM_=uo{jX#qbCb6>_1WN_HKjrpZr3;5|fxIwu>%xgiN64r2yBV+B1PrGTJ-Ha^ zLh`QgBPK6iWh_eC^Vcp5kmA3a|2uXHcR%o%6gmDn7!GYL6vx?{Yse#z-SlS+oGs?$DW*3#&xv3KGw&Uxu+A^yBuQZ0}7zlRb& z;u)PlRHE_tB%KLVJA`n)-zb?T6*=c4-qH?I{S7(_-e2Mt2VqY^?y^@-8H`901V?v) zg{(TB&UWUK0A``EQyhT5`wVHsVjH;lJO=(A-0VugHFNN)=a=JkH1Iy2m8O*|zGW`% zL^29wKEcnmKQYYk7ms&MB~WQa$Zt5@`TPqFz*DI@&gH5USaB~ZGyUhKGXEyoQqSBx z|Ikvm(t4Y*AW0ldnsWy_F2Jw_fUkQ*9 zTiltoog`Bq{MeF6Qm`{66J(dba9ZNr@6t$zpK7p=qD4os#h!*fiGj&JfIJ@k*qFXTC=#W#dpdCP9 zEQaKDW$}PZD9HIIhA?Q-YktJjLezhK%ku7YcKtCxXd<{x9q+GhrDNcM} zoKk46V!{W$M!#e9K(?&*$E{AUu4HN~T9_DAGS;?~Vfl`$s3G-|5MP^gOVdqEvm*(DZ3zAz9ezI5O33cMe;1Q4o8+;b4OF z`nY)E&SN{_WY{zRQ(21#f5nOabqupNbBm{^x#<#!QP;T`#w0~q(4kL4+tGEw!h9*8 zroC0tLIjmbPd7-?0%~!rr?A~RHiz$QF9_2SX;%gmYzsFbKyve%`G$$I2ibm#9Jro& z{OubC%X!=Qo|@|`M$R8(GZRI*R{Je2n(2Ioo}nh!Zx_B0jO0h{lTH+;2OmnKU+3e( z#lKTvN#$AHnVLH{8P~hGU9tEZmPA>lAV2@kqmO!9gELbK^QhC&&EFyj?E*3ctI+Yb zofKKp^E`;fsjT3dLVuuLa8uw_YoBaedE)!Zlp=C-1G!*E%h$gaxG_xj!%>D?4haX` z&SIh&`}kt_X@7#pz_2J}9|!W+u_*;QsOmi+MCD)N8M2ge>Flz-D_*3&V@&3>wzEM# zK*mb_qp)4S;HoJM0cP%Xx5NXXJ!^X1U|7jP2%e+Y_pDeFl5|uN!(aL1!zu0=pPg8< ze6+(KFV-*x1tvf-#);9evRwG&srT43V)xr87<_&kAqX3MX!W%Bws!~j2X0NJdv_Wj zn(l`+F^bRc4DoFQIu$#-GZ^8M?VJFi>-OxL(ZDeC2aJ*|1SAFk0DuMjC!g{Y z&-gbp0R-S#2MU1pb=S_w)z-+EIoOvu+?N^b3wWC^GlUqyf3h&c`(z=$Kz@Px0{sQ% z3+xxTFYsRwz94=<`hxrgQvVg76f@T=Aa{sb zi7G+rAl5GEhyKjc?-EJ8)obe_=1@C=L>B+ofe{w@r!SchQH%MHIM`eqoN8h4W0=v| z02qdcwB-@%1x$52isU~dVE-V@B|udazNYs!%_okUg{wMZ%9LJX&Ue)xS|_j&25OZ* zhl;k_HT8}ff}1++s_Og?okx=Zwax9(-YVWOD{d1A`yzj*Z8-e!xw2tlUK38(&8sZ4 zzV~~xGGKXnFl~>sv(v%7qAVmS!<2OtKdZ6JFr+*EY3U70&4}CQ+cR^bJyqgC42I?o z_yXle)4&UUTfq@|7{q@bS32z(x%B@(oi-hKM$Yl`|F$B&76kz8KLtR6`Tt!X4sNcF zZmvw8cD72gVBi42|F@trK)vd@N5lXU1QG!1e~0}~QIn2#PEJhxpE=_7!70{h004Sc z!2i+7BVF;U(|?0>Iw!{m3rdgupLX4U+S=*k=j7ys|9#V6QK5d#jsO2l@V|qUIwViG tKO`qer@kPk15Z9m&$}SUPq)7yX97!3&%Gd5{kmLVkn=+)pOSwK{9jAy&rSdU diff --git a/functions/getIntermediateU_xyRK3_dst.m b/functions/getIntermediateU_xyRK3_dst.m new file mode 100644 index 0000000..3cde381 --- /dev/null +++ b/functions/getIntermediateU_xyRK3_dst.m @@ -0,0 +1,30 @@ +% Copyright 2020 The MathWorks, Inc. +function xu = getIntermediateU_xyRK3_dst(u, b, dt, Re, nx, ny, dx, dy, id) + +kx = [1:nx-1]'; +ax = pi*kx/(nx); +mwx = 2*(cos(ax)-1)/dx^2;% DST-I + +ky = [1:ny]'; +ay = pi*ky/(ny); +mwy = 2*(cos(ay)-1)/dy^2; % DST-II + +mw = mwx+mwy'; % Modified Wavenumber + +% 各 column 毎に変換されるので、 +% column 方向(x 方向に)dst1, 転置して +% y 方向に dst2、転置し返して元に戻す処理 +rhshat = mydst2(mydst1(b)')'; + +switch id + case 1 + uhat = rhshat./(1-(4/15)*dt/(Re)*mw); + case 2 + uhat = rhshat./(1-(1/15)*dt/(Re)*mw); + case 3 + uhat = rhshat./(1-(1/6)*dt/(Re)*mw); +end + +xu = mydst1(mydst3(uhat')'); + +end \ No newline at end of file diff --git a/functions/getIntermediateV_xyRK3_dst.m b/functions/getIntermediateV_xyRK3_dst.m new file mode 100644 index 0000000..2fbbcfe --- /dev/null +++ b/functions/getIntermediateV_xyRK3_dst.m @@ -0,0 +1,32 @@ +% Copyright 2020 The MathWorks, Inc. +function xv = getIntermediateV_xyRK3_dst(v, b, dt, Re, nx, ny, dx, dy,id) + + +kx = [1:nx]'; +ax = pi*kx/(nx); +mwx = 2*(cos(ax)-1)/dx^2;% DST-I + +ky = [1:ny-1]'; +ay = pi*ky/(ny); +mwy = 2*(cos(ay)-1)/dy^2; % DST-II + +mw = mwx+mwy'; % Modified Wavenumber + +% 各 column 毎に変換されるので、 +% column 方向(x 方向に)dst2, 転置して +% y 方向に dst1、転置し返して元に戻す処理 +rhshat = mydst1(mydst2(b)')'; + +switch id + case 1 + vhat = rhshat./(1-(4/15)*dt/(Re)*mw); + case 2 + vhat = rhshat./(1-(1/15)*dt/(Re)*mw); + case 3 + vhat = rhshat./(1-(1/6)*dt/(Re)*mw); +end + +xv = mydst3(mydst1(vhat')'); + +end + diff --git a/functions/updateVelocityField_RK3.m b/functions/updateVelocityField_RK3.m index b7afc4f..6b96d19 100644 --- a/functions/updateVelocityField_RK3.m +++ b/functions/updateVelocityField_RK3.m @@ -81,10 +81,12 @@ % Implicit treatment for xy direction b = u(2:end-1,2:end-1) - dt*(gamma(id)*Nu + zeta(id)*Nu_old - alpha(id)/(Re)*(Lux+Luy+Lubc)); -xu = getIntermediateU_xyRK3(u, b, dt, Re, nx, ny, dx, dy, id); +% xu = getIntermediateU_xyRK3(u, b, dt, Re, nx, ny, dx, dy, id); +xu = getIntermediateU_xyRK3_dst(u, b, dt, Re, nx, ny, dx, dy, id); b = v(2:end-1,2:end-1) - dt*(gamma(id)*Nv + zeta(id)*Nv_old - alpha(id)/(Re)*(Lvx+Lvy+Lvbc)); -xv = getIntermediateV_xyRK3(v, b, dt, Re, nx, ny, dx, dy, id); +% xv = getIntermediateV_xyRK3(v, b, dt, Re, nx, ny, dx, dy, id); +xv = getIntermediateV_xyRK3_dst(v, b, dt, Re, nx, ny, dx, dy, id); u(2:end-1,2:end-1) = xu; v(2:end-1,2:end-1) = xv;