A short review of how to change the properties of lattice node

This post was written in English because my workstation(contains source codes) was not installed with any Chinese input method. In the next few days I’ll install one.

Cell Location

By the code below you can reach to the cell of the location you want as “cell”. Then it will be easier to make other operations.

Cell<T,Descriptor>& cell  = lattice.get(ix,iy,iz);
For parallelism

If you would like to run the code in parallel, the below operation is a necessity.

for (plint iX=domain.x0; iX<=domain.x1; ++iX) {
for (plint iY=domain.y0; iY<=domain.y1; ++iY) {
for (plint iZ=domain.z0; iZ<=domain.z1; ++iZ) {
Cell<T,Descriptor>& cell = lattice.get(iX,iY,iZ);
......
}}}

If you would like to get the global coordinates, you can refer to the code below

Dot3D absoluteOffset = lattice.getLocation();
for (plint iX=domain.x0; iX<=domain.x1; ++iX)
{
plint absoluteX = absoluteOffset.x + iX;
for (plint iY=domain.y0; iY<=domain.y1; ++iY)
{
plint absoluteY = absoluteOffset.y + iY;
for (plint iZ=domain.z0; iZ<=domain.z1; ++iZ)
{
plint absoluteZ = absoluteOffset.z + iZ;
Cell<T, Descriptor>& cell = lattice.get(iX,iY,iZ);

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

}
}
}

Code source: https://palabos-forum.unige.ch/t/data-processor-for-non-local-operation/1075/2

iniCellAtEquilibrium(cell, density, velocity)
T density = cell.computeDensity(); 
Array<T, 3> velocity;
cell.computeVelocity(velocity); 
		  cell.defineDensity(density);
		  cell.defineVelocity(velocity);
          iniCellAtEquilibrium(cell, density, velocity);	

By using function inicellAtEquilibrium, you will need density and velocity to assign the Feq to that cell. By cell.computeDensity() you will get the density and by cell.computeVelocity(velocity) you will calculate that cell’s velocity and assign it to the variable you defined, for example, here is the “velocity”.

So you want to change the properties of the fluid node during simulation, you can add codes to modify the density and velocity before defineDensity and defineVelocity.

computeEquilibria(fEq, rhoBar, j, jSqr)
template<typename T, template<typename U> class Descriptor> void Dynamics<T,Descriptor>::computeEquilibria (
        Array<T,Descriptor<T>::q>& fEq,  T rhoBar, Array<T,Descriptor<T>::d> const& j, T jSqr, T thetaBar ) const {
    for (int iPop=0; iPop<Descriptor<T>::q; ++iPop) {
        fEq[iPop] = computeEquilibrium(iPop, rhoBar, j, jSqr, thetaBar);
    } }

Before I introduce the computeEquilibria, the above code of function should be demonstrated. It shows that the variable fEq is calculated by function computeEquilibriumwith inputs iPop, rhoBar, j, and jSqr, where iPop depends on the lattice structure you chose. For example if you use D3Q19, then the ::q term here will be 19.

Then we look at the main codes to see how this way can perform.

Cell<T,Descriptor>& cell = lattice.get(iX,iY,iZ);
T rhoBar; Array<T,3> j; Array<T,D::q> fEq;
cell.getDynamics().computeRhoBarJ(cell, rhoBar, j);
T jSqr = normSqr(j);
cell.getDynamics().computeEquilibria(fEq, rhoBar, j, jSqr);

The idea is simple, firstly calculate the rho, j, and jsqr, and then we will gather everything we need to get the Feq. The computeEquilibria only does the favor of looping for us.

If we want to change the properties of that cell, we can modify rhoBar, j to assign the new Feq to fEq and implement the code below.

for (plint iPop=0; iPop<D::q; ++iPop) {
		cell[iPop] = fEq[iPop];
}
computeEquilibrium(iPop, rhoBar, j, jSqr)

This way is similar to the method above, I think I don’t need to make further explanation.

for(plint iPop=0;iPop<D::q;++iPop){
cell[iPop] = cell.computeEquilibrium(iPop, rhoBar, j, jSqr);
}
To calculate the phase change energy from bcc to fcc using LAMMPS, you can use the lattice orientation relationship method. This method involves creating a bcc and fcc crystal in LAMMPS, aligning their lattice orientations, and then calculating the energy difference between the two structures. Here is an example LAMMPS input script that demonstrates this method: ``` # Create bcc crystal units metal atom_style atomic lattice bcc 2.855 region box block 0 10 0 10 0 10 create_box 1 box create_atoms 1 box mass 1 55.85 # Relax bcc structure pair_style eam pair_coeff * * Fe_u3.eam neighbor 2.0 bin neigh_modify every 20 delay 0 check no fix 1 all box/relax iso 0.0 vmax 0.001 thermo 10 thermo_style custom step pe lx ly lz press pxx pyy pzz minimize 1e-8 1e-8 1000 10000 # Create fcc crystal clear units metal atom_style atomic lattice fcc 3.615 region box block 0 10 0 10 0 10 create_box 1 box create_atoms 1 box mass 1 55.85 # Align lattice orientations clear units metal atom_style atomic lattice fcc 3.615 orient x 1 0 0 orient y 0 1 0 orient z 0 0 1 region box block 0 10 0 10 0 10 create_box 1 box create_atoms 1 box mass 1 55.85 group fcc type 1 group bcc type 2 fix 1 fcc setforce 0.0 0.0 0.0 fix 2 bcc setforce 0.0 0.0 0.0 minimize 1e-8 1e-8 1000 10000 variable dx equal lx/10.0 variable dy equal ly/10.0 variable dz equal lz/10.0 displace_atoms bcc move ${dx} 0.0 0.0 # Calculate energy difference compute bccpe bcc pe compute fccpe fcc pe variable dE equal c_fccpe-c_bccpe print "Energy difference: ${dE} eV/atom" ``` This script first creates a bcc crystal and relaxes its structure using the EAM potential for iron. It then creates an fcc crystal and aligns its lattice orientation with the bcc crystal. The two structures are then minimized to their equilibrium positions. Finally, the energy difference between the fcc and bcc structures is calculated using the compute command in LAMMPS. Note that you will need to replace the EAM potential file `Fe_u3.eam` with the appropriate potential for your system. I hope this helps!
评论 4
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值