关于insert和periodic的问题

前言 

​​​​​​​各种宏的定义及作用_imastrid的博客-CSDN博客

在我前面写的这一篇里面，我提到了mppic代码中configure.h文件的设置宏的定义和作用。

最近才好像有点搞懂了代码中是怎么插入parcel或者循环parcel的，所以这里想整理一下突然明白的部分。

定义及意义

不厌其烦地再提起一次定义及意义。

__INSERT_PARCEL_

它是用来定义采用插入的方法来加入颗粒。

代码中具体实现方法是：每一个步长在固定的位置或者在通过一些判断和随机数确定的位置（比如入口处或者分布板下面或者侧边入口）加入定量或者随机数的parcel数（按照一些工业实际量计算）。

针对MIP反应器实际工业量的insert量可以通过这个matlab代码来计算：matlab计算按步长加入的颗粒量_imastrid的博客-CSDN博客

__PERIODIC__

它是用来实现循环颗粒的作用的。比如颗粒从出口出，那就将它回收到入口。

后面也可以自己再改，比如将入口部分边界条件直接等于出口的边界条件。

代码中具体实现方法：判断颗粒是否在domain里面（ps[id].in=1,0,-1,-2,-3，1表示在，其它情况不清楚），如果不在，就将它回收到reinsert.ini（init/文件夹中）里面定义的区域里，速度重新定义为0或者任何其他数。

代码实现

__INSERT_PARCEL_

对于Ps=010的串行代码：ps[id].in只有1，0，-1

出现位置：configure.h（定义宏），mppic_main.cpp→insert.cpp，create_parcel.cpp

1. create_parcel.cpp

里面的give_property()函数，定义完ps（parcelStay）的颗粒数，直径以及体积，对于__INSERT_PARCEL__初始化给ps[id].in=-1，否则没有__INSERT_PARCEL__则为1。

2. mppic_main.cpp→insert.cpp

在这里面需要手动定义了插入的位置和插入数量，还有选择插入的方式（有cuboid，column，cirque），最后记得update_index

对于TEMP_MPI3的并行代码：ps[id].in有1，-1，-2，-3

出现位置：configure.h（定义宏），mppic_main.cpp→insert.cpp，create_parcel.cpp（但被注释掉了，因此跟这个宏没关系其实，但放在这里与上面做对比）

1. create_parcel.cpp

里面的give_property()函数，定义完ps（parcelStay）的颗粒数，直径以及体积，ps[id].in=-3。

2. mppic_main.cpp→insert.cpp

在这里面需要手动定义了插入的位置和插入数量，还有选择插入的方式（有cuboid，column，cirque），最后记得update_index

__PERIODIC__

对于Ps=010的串行代码：ps[id].in只有1，0，-1

出现位置：configure.h（定义宏），mppic_main.cpp→read_reins.cpp，walls.cpp，periodicReinsert.H，

1. mppic_main.cpp→read_reins.cpp

读取reinsert.ini里面的数据。

2. walls.cpp

如果parcel跑出了outlet边界（用otri判断），如果是__PERIODIC__则定义ps[id].in=0，否则没有__PERIODIC__则为-1

3. periodicReinsert.H

对于__PERIODIC__，判断ps[id].in=0的颗粒是否outside simulation domain，如果是就根据reinsert.ini定义的区域和速度重新赋值，重新赋值后就是1。没有__PERIODIC__就是只有1和-1。

要么只有1和-1，要么有1，0，-1，这个我觉得有点奇怪-1没回收，这不就是跑了啊？

对于TEMP_MPI3的并行代码：ps[id].in有1，-1，-2，-3

出现位置：configure.h（定义宏），mppic_main.cpp→read_reins.cpp，walls.cpp，periodicReinsert.H，update.cpp

1. mppic_main.cpp→read_reins.cpp 

读取reinsert.ini里面的数据。

2. walls.cpp

如果parcel跑出了outlet边界（用otri判断），如果是__PERIODIC__则定义ps[id].in=-2，否则-1

3. periodicReinsert.H

（这一段被掉屏蔽了，主函数也没调用这个头文件，因此内容仍与串行代码一致，没改成-2，我感觉如果想将它跑成串行则应该是要将0改成-2的）对于__PERIODIC__，判断ps[id].in=0的颗粒是否outside simulation domain，如果是就根据reinsert.ini定义的区域和速度重新赋值，重新赋值后就是1。

4. update.cpp

update.cpp中的mpi_data_trans_parcels()函数中，对于__PERIODIC__，添加了最上面的processor和最底下的processor的信息交换（ParcelTran）的代码（这也说明了，这个代码并行只能采用垂直切割的方法，其他方法就传递不了准确信息了），例子：

#ifdef __PERIODIC__
        MPI_Barrier(MPI_COMM_WORLD);
        if (proc.rank == proc.size - 1) {
            MPI_Send(&nptf, 1, MPI_INT, 0, dim * 2, MPI_COMM_WORLD);
        }
        if (proc.rank == 0) {
            MPI_Recv(&npfb, 1, MPI_INT, proc.size - 1, dim * 2, MPI_COMM_WORLD,
                    &proc.status);
        }

        MPI_Barrier(MPI_COMM_WORLD);
#endif

mpi_data_trans_parcels()函数

1. 首先数一数有多少不是ps[id].in=-3和-2，并且位置不在本区域内部的。
2. 然后将nptf数据地址发送给npfb，或者nptb发给npff。指针都指向相同地址时，数据就一样了。
3. 然后定义一维的nptf或者nptb个ParcelTran结构ptf,pfb,ptb,pff。 
4. 如果ps[id].in=-3和-2，位置不在本区域内部的，赋值pt数据给定义的ParcelTran结构，另外定义ps[id].in=-2。
5. 然后将ptf传递信息给pfb，将ptb传递给pff，对于PERIODIC，则要再更新一下出口到进口循环回来的颗粒纵坐标。
6. 传递完以后，将所有信息传回给pt，然后ps[id].in!=1的都定义成1，然后npeach++，算一下多少在区域内。
7. 同时用npff,npfb,nptf.nptb更新prmt.npp，清空一维结构。

ps[id].in不同数值的意义

串行

1：初始化无insert时所有parcel为1；在区域内

0：定义了periodic的不在区域内

-1：初始化时如果是insert则所有parcel初始为-1；不在区域内

串行由于它没用insert，因此它的逻辑不适用。初始化全是-1，后面都运行不起来了。

但是periodic的逻辑可以：

1. create_parcel初始化为1，计算
2. 然后到periodicReinsert的时候，判断再分类1和-1
3. 然后到mppic_update_one_step的时候，只对1的颗粒计算，update_1st和update_2nd中的stl_wall_collide_1st(id)和stl_wall_collide_2nd(id)，如果颗粒穿出outlet则为0。
4. 然后这些0在periodicReinsert的时候都会被回收循环。
5. 对于-1好像就没有什么操作，判断完以后它们就没了等于。

并行

1：初始化1，在区域内

-1：出现在很多被屏蔽的walls.cpp中的语句，还出现在periodicReinsert.H的屏蔽语句中；对于未屏蔽的语句，反正就是不在区域内时无特殊定义就为-1。

-2：walls.cpp屏蔽语句里面，如果定义periodic就是-2；mpi_data_trans_for_parcels里面在将pt值传到缓冲区ptf和ptb时，设为-2，然后将所有pfb和pff传给pt时设回1。

-3：初始化give_property时-3

其中，感觉-2的颗粒是可以回收变成1的，-1的颗粒都不要了，-3的颗粒也是不要的，我甚至不知道-3的意义，可能可以用来计算相对初始少了多少parcel吧。

并行代码总定义了macros.h里面有针对-3和-2要continue的宏。

很多句子都被屏蔽掉了，只针对TEMP_MPI3案例的话，有1，-1，-2，-3出现在未屏蔽语句中。

由于是insert，则其逻辑为：

1. create_parcel读init/vtk文件，give_property初始化-3；然后读result/vtk文件，对nptmp个颗粒赋值，还要判断速度不超过某些定值，最后read_vtk2更新成1。剩下的-3就是init文件中的不需要的。（虽然我不知道这个定义-3的意义是啥）
2. 然后到mppic_update_one_step的时候，首先是insert_parcel()，所有insert的都是1，接着只对1的颗粒计算，对于update里面有mpi_data_trans_for_parcels()函数，将pt值传到缓冲区ptf和ptb时，设为-2，然后将所有pfb和pff传给pt时设回1；下一步check_all_parcels的时候会将所有1再次分类1和-1，但没有再对穿出outlet的颗粒赋值了。

使用方法

所有代码具体选项应该是我们针对不同计算案例而手动去修改的，因此了解代码的位置和修改方法很重要。上面已经完成了这些总结。下面简单介绍一下，怎么使用：

• 注意宏定义，在configure.H文件中，需要哪一部分就修改哪一部分代码，屏蔽不要的代码
• 跑案例前记得检查逻辑，关于初始化ps[id].in的值的逻辑要检查
• 对于反应器模拟，串并行的区别要注意，特别是它们对periodic的处理是分开的，不在同一个地方
• 如果是定量加入催化剂，则采用insert，如果是全床颗粒量不变，则采用periodic。在采用insert的时候算出一个稳定状态时，则可以根据这个算例算一下全床平均孔隙率，之后可以采用periodic快速达到稳定状态。

....................................................................................................................................................

20211220更

突然明白了periodicReinsert.H的意义，它可以拿来收集不在区域内跑到边界外的，也可以一次性收集从出口出来的。总之只要一个案例中所有函数一起做到把跑出去的（不明原因跑到边界外的，从出口出去的，而且要注意串并行的区别）都收集了就行了。