MFiX利用ReactionRates数组把颗粒变量存储到网格上

思路

受到这篇Jay的文章的启发，修改了一下
https://www.cnblogs.com/Jay-CFD/p/10786299.html
借鉴原文的方法二

但是这个方法修改了vtk_out.f这个文件，由于19版本和20版本的MFiX在这个文件上改动很大，所以在20上无法运行成功。

于是做了如下修改：
我发现原博主更改vtk_out.f这个文件主要其实就是用来call des_granular_temperature的，所以何不在别的地方call des_granular_temperature呢？

由此进一步搜索具有call des_granular_temperature的代码，发现
/legacy_tests/dem-tests/homogeneous-cooling-system/usr0_des.f这个测试算例调用了des_granular_temperature。

于是复制到当前文件夹下。并在GUI中勾选call usr routine。编译后运行，发现想法验证成功！

下面是具体操作流程：

操作流程

原文的例子是要输出颗粒的Y速度到网格内。使用ReactionRates这个数组。
注：所有修改的源文件先复制到自己的文件夹下

1. 修改des_granular_temperature.f

这部分和Jay的原文一样
一共五个地方修改

!##################CL_200917##################
   use rxns,only: reactionrates,nrr
!#############################################

!#####################CL_200917###############
      DOUBLE PRECISION :: v_s_tmp(DIMENSION_3, DIMENSION_M)
!#############################################

!#####################CL_200917###############
      v_s_tmp(:,:)=ZERO
!#############################################

!#####################CL_200917###############
      v_s_tmp(ijk,M)=v_s_tmp(ijk,M) + DES_VEL_NEW(LL,2)
!#############################################

!#####################CL_200917###############
                    reactionrates(ijk,1)=v_s_tmp(ijk,M)/&
                    DBLE(NP_PHASE(IJK,M))
!#############################################

2. 将homogeneous-cooling-system/usr0_des.f复制到当前文件夹

具体位置在
/legacy_tests/dem-tests/homogeneous-cooling-system/usr0_des.f

!vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvv!
!                                                                      !
!  Module name: URS0_DES                                               !
!                                                                      !
!  Purpose: This routine is called before the discrete phase time loop !
!  and is user-definable. The user may insert code in this routine or  !
!  call appropriate user defined subroutines.                          !
!                                                                      !
!  This routine is not called from a loop, hence all indicies are      !
!  undefined.                                                          !
!                                                                      !
!  Author: J.Musser                                   Date: 06-Nov-12  !
!                                                                      !
!  Comments:                                                           !
!                                                                      !
!^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^!
      SUBROUTINE USR0_DES

      Use des_rxns
      Use des_thermo
      Use discretelement
      Use run
      Use usr

      IMPLICIT NONE

      CALL WRITE_GRAN_TEMP

      RETURN
      END SUBROUTINE USR0_DES


!vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvv!
!                                                                      !
!  Module name: URS0_DES                                               !
!                                                                      !
!  Purpose: This routine is called before the discrete phase time loop !
!  and is user-definable. The user may insert code in this routine or  !
!  call appropriate user defined subroutines.                          !
!                                                                      !
!  This routine is not called from a loop, hence all indicies are      !
!  undefined.                                                          !
!                                                                      !
!  Author: J.Musser                                   Date: 06-Nov-12  !
!                                                                      !
!  Comments:                                                           !
!                                                                      !
!^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^!
      SUBROUTINE WRITE_GRAN_TEMP

      use des_granular_temperature_mod, only: des_granular_temperature
      use discretelement, only: DES_MMAX
      use discretelement, only: DES_VEL_AVG
      use discretelement, only: GLOBAL_GRAN_TEMP
      use discretelement, only: OVERLAP_MAX
      use discretelement, only: S_TIME
      use discretelement, only: des_ke
      use geometry, only: NO_K
      use run, only: RUN_NAME

      IMPLICIT NONE

! Local variables
!---------------------------------------------------------------------//
! file name
      CHARACTER(LEN=64) :: FNAME
! logical used for testing is the data file already exists
      LOGICAL :: F_EXISTS
! file unit for heat transfer data
      INTEGER, PARAMETER :: GT_UNIT = 2030

! Variables used for calculation in output file
      DOUBLE PRECISION :: GRAN_TEMP, AVG_VEL
! Working dimension
      INTEGER :: wDIMN

! Set the working dimension for the problem.
      wDIMN = merge(2, 3, NO_K)

      CALL DES_GRANULAR_TEMPERATURE

      GRAN_TEMP = (1.0/wDIMN)*SUM( GLOBAL_GRAN_TEMP(1:wDIMN) )
      AVG_VEL = (1.0/wDIMN)*SUM( DES_VEL_AVG(1:wDIMN) )


      FNAME = 'POST_GT.dat'
      INQUIRE(FILE=FNAME,EXIST=F_EXISTS)

      IF(.NOT.F_EXISTS) THEN
         OPEN(UNIT=GT_UNIT,FILE=FNAME,STATUS='NEW')
         WRITE(GT_UNIT,"(4(2X,A))") 'SOLIDS-TIME', 'GRAN-ENERGY', &
            'KNTC-ENERGY', 'AVG-VEL'
      ELSE
         OPEN(UNIT=GT_UNIT,FILE=FNAME,&
            POSITION="APPEND",STATUS='OLD')
      ENDIF

      WRITE(GT_UNIT,1200) S_TIME, GRAN_TEMP, DES_KE, AVG_VEL

 1200 FORMAT(4(2x,g11.5))

      CLOSE(GT_UNIT)

      RETURN
      END SUBROUTINE WRITE_GRAN_TEMP

3 修改算例.mfx文件

算例采用的还是fluid_bed_dem_2d这个官方算例
new project里面直接就能导入了

但是要 作如下几点修改(除了第一点以外其他都和原文一致）

1. 勾选enable user defined subroutine(比Jay原文新增)
2. 给定ReactionRates array size
   
3. 在output中新增输出网格的vtk，并指定ReactionRates变量名称
   

算例结果验证

可以看到输出结果了。如图