发现粘性力功率计算的实现和公式(4.59)有出入 #11
Closed
Conversation
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
|
@ZZP-DMU，感谢提出意见和修改，我将尽可能在这周内抽空探查所述问题。 |
zoziha
approved these changes
Dec 5, 2023
| ! tdsdt(j) = tdsdt(j) + txx(i)*txx(j) + 2._rk*txy(i)*txy(j) + 2._rk*txz(i)*txz(j) + & | ||
| ! tyy(i)*tyy(j) + 2._rk*tyz(i)*tyz(j) + tzz(i)*tzz(j) | ||
| ! end if | ||
| ! end do |
There was a problem hiding this comment.
但由于这里只会影响能量守恒公式的计算结果（内能），不会影响程序的加速度求解，不知道粘性熵的具体公式，建议保持原状，保留你的注释。
There was a problem hiding this comment.
tdsdt似乎是指粘性力功率，粘性力的公式书中在能量方程中分离出来了。不知道是否是指这个
src/single_step/internal_force.f90
Outdated
| ! tyy(i)*tyy(j) + 2._rk*tyz(i)*tyz(j) + tzz(i)*tzz(j) | ||
| ! end if | ||
| ! end do | ||
| ! tdsdt(i) = 0.5_rk*eta(i)/rho(i)*tdsdt(i) |
There was a problem hiding this comment.
这段应该是：
do i = 1, ntotal
tdsdt(i) = 0.5_rk*eta(i)/rho(i)*tdsdt(i)
end do| @@ -53,6 +53,7 @@ subroutine kernel(r, dx, hsml, w, dwdx) | |||
| w = factor*1._rk/6._rk*(2._rk - q)**3 | |||
| do d = 1, dim | |||
| dwdx(d) = -factor*1._rk/6._rk*3.*(2._rk - q)**2/hsml*(dx(d)/r) | |||
| !dwda(d) =-factor*(2._rk - q)**2/2*hsml*(dx(d)/r) | |||
| @@ -82,6 +82,7 @@ subroutine direct_find(itimestep, ntotal, hsml, x, niac, pair_i, pair_j, w, dwdx | |||
| !> kernel and derivations of kernel | |||
|
|
|||
| call kernel(r, dxiac, mhsml, w(niac), tdwdx) | |||
| !kernel(r, dx, hsml, w, dwdx) | |||
| do d = 1, dim | |||
| @@ -156,32 +156,62 @@ subroutine int_force(itimestep, dt, ntotal, hsml, mass, vx, niac, rho, eta, pair | |||
| hvcc = 0._rk | |||
| do d = 1, dim | |||
| hvcc = hvcc + dvx(d)*dwdx(d, k) | |||
| !看不懂hvcc&vcc到底是在计算什么东西 | |||
There was a problem hiding this comment.
在/src/sing_step/art_heat.f90中vcc&hvcc都出现了，应该和这里类似，表示速度的散度。但这里确实没有使用vcc&hvcc
src/single_step/density.f90
Outdated
| do i = 1, ntotal | ||
| rho(i) = rho(i)/wi(i) | ||
| end do | ||
| end if | ||
| !？归一化，不太理解？ |
src/single_step/external_force.f90
Outdated
| @@ -52,6 +52,14 @@ subroutine ext_force(ntotal, mass, x, niac, pair_i, pair_j, itype, hsml, dvxdt) | |||
| do i = 1, ntotal | |||
| dvxdt(dim, i) = -9.8_rk | |||
There was a problem hiding this comment.
这里可能写成：
dvxdt(dim, i) = dvxdt(dim, i) - 9.8_rk但意义差不多，问题不大。
There was a problem hiding this comment.
是的，dvxdt仅是输出的参数。
我之前的理解是：如果从外部传入的话，dvxdt可能会有非零初始值。
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment
Add this suggestion to a batch that can be applied as a single commit.
This suggestion is invalid because no changes were made to the code.
Suggestions cannot be applied while the pull request is closed.
Suggestions cannot be applied while viewing a subset of changes.
Only one suggestion per line can be applied in a batch.
Add this suggestion to a batch that can be applied as a single commit.
Applying suggestions on deleted lines is not supported.
You must change the existing code in this line in order to create a valid suggestion.
Outdated suggestions cannot be applied.
This suggestion has been applied or marked resolved.
Suggestions cannot be applied from pending reviews.
Suggestions cannot be applied on multi-line comments.
Suggestions cannot be applied while the pull request is queued to merge.
Suggestion cannot be applied right now. Please check back later.
公式(4.59)中计算粘性力功率的项是：
$\frac{\mu_i}{2 \rho_i}\varepsilon_i^{\alpha \beta} \varepsilon_j^{\alpha \beta} $
但是在/src/single_step/internal_force.f90中的粘性力功率实现代码是:
我认为应该修改为如下形式