Commit e2c84bd4 authored by user's avatar user
Browse files

Changes to be committed: 

	modified:   coefs.h
	modified:   core.cpp
	modified:   core.h
	modified:   lbm.cpp
	modified:   main.cpp
	modified:   taylor.cpp
	modified:   taylor.h
	new file:   tex/buildall
	new file:   tex/d1q3/d1q3_C.tex
	new file:   tex/d1q3/d1q3_defs.tex
	new file:   tex/d1q3/d1q3_matM.tex
	new file:   tex/d1q3/d1q3_matM_per_els.tex
	new file:   tex/d1q3/d1q3m0_feq.tex
	new file:   tex/d1q3/d1q3m0_matA.tex
	new file:   tex/d1q3/d1q3m0_matA_per_els.tex
	new file:   tex/d1q3/d1q3m0_matB.tex
	new file:   tex/d1q3/d1q3m0_matB_per_els.tex
	new file:   tex/d1q3/d1q3m0_matM.tex
	new file:   tex/d1q3/d1q3m0_matM_per_els.tex
	new file:   tex/d1q3/d1q3m1_feq.tex
	new file:   tex/d1q3/d1q3m1_matA.tex
	new file:   tex/d1q3/d1q3m1_matA_per_els.tex
	new file:   tex/d1q3/d1q3m1_matB.tex
	new file:   tex/d1q3/d1q3m1_matB_per_els.tex
	new file:   tex/d1q3/d1q3m1_matM.tex
	new file:   tex/d1q3/d1q3m1_matM_per_els.tex
	new file:   tex/d1q3/d1q3m5_feq.tex
	new file:   tex/d1q3/d1q3m5_matA.tex
	new file:   tex/d1q3/d1q3m5_matA_per_els.tex
	new file:   tex/d1q3/d1q3m5_matB.tex
	new file:   tex/d1q3/d1q3m5_matB_per_els.tex
	new file:   tex/d1q3/d1q3m5_matK.tex
	new file:   tex/d1q3/d1q3m5_matK_per_els.tex
	new file:   tex/d1q3/d1q3m5_matM.tex
	new file:   tex/d1q3/d1q3m5_matM_per_els.tex
	new file:   tex/d1q3/supp_d1q3_defs.tex
	new file:   tex/d1q3/supp_d1q3m0_defs.tex
	new file:   tex/d1q3/supp_d1q3m1_defs.tex
	new file:   tex/d1q3/supp_d1q3m5_defs.tex
	new file:   tex/d2q5/d2q5_C.tex
	new file:   tex/d2q5/d2q5_defs.tex
	new file:   tex/d2q5/d2q5_matM.tex
	new file:   tex/d2q5/d2q5_matM_per_els.tex
	new file:   tex/d2q5/d2q5m0_feq.tex
	new file:   tex/d2q5/d2q5m0_matA.tex
	new file:   tex/d2q5/d2q5m0_matA_per_els.tex
	new file:   tex/d2q5/d2q5m0_matB.tex
	new file:   tex/d2q5/d2q5m0_matB_per_els.tex
	new file:   tex/d2q5/d2q5m0_matM.tex
	new file:   tex/d2q5/d2q5m0_matM_per_els.tex
	new file:   tex/d2q5/d2q5m1_feq.tex
	new file:   tex/d2q5/d2q5m1_matA.tex
	new file:   tex/d2q5/d2q5m1_matA_per_els.tex
	new file:   tex/d2q5/d2q5m1_matB.tex
	new file:   tex/d2q5/d2q5m1_matB_per_els.tex
	new file:   tex/d2q5/d2q5m1_matM.tex
	new file:   tex/d2q5/d2q5m1_matM_per_els.tex
	new file:   tex/d2q5/d2q5m2_feq.tex
	new file:   tex/d2q5/d2q5m2_matA.tex
	new file:   tex/d2q5/d2q5m2_matA_per_els.tex
	new file:   tex/d2q5/d2q5m2_matB.tex
	new file:   tex/d2q5/d2q5m2_matB_per_els.tex
	new file:   tex/d2q5/d2q5m2_matM.tex
	new file:   tex/d2q5/d2q5m2_matM_per_els.tex
	new file:   tex/d2q5/d2q5m3_feq.tex
	new file:   tex/d2q5/d2q5m3_matA.tex
	new file:   tex/d2q5/d2q5m3_matA_per_els.tex
	new file:   tex/d2q5/d2q5m3_matB.tex
	new file:   tex/d2q5/d2q5m3_matB_per_els.tex
	new file:   tex/d2q5/d2q5m3_matM.tex
	new file:   tex/d2q5/d2q5m3_matM_per_els.tex
	new file:   tex/d2q5/d2q5m4_feq.tex
	new file:   tex/d2q5/d2q5m4_matA.tex
	new file:   tex/d2q5/d2q5m4_matA_per_els.tex
	new file:   tex/d2q5/d2q5m4_matB.tex
	new file:   tex/d2q5/d2q5m4_matB_per_els.tex
	new file:   tex/d2q5/d2q5m4_matM.tex
	new file:   tex/d2q5/d2q5m4_matM_per_els.tex
	new file:   tex/d2q5/d2q5m5_feq.tex
	new file:   tex/d2q5/d2q5m5_matA.tex
	new file:   tex/d2q5/d2q5m5_matA_per_els.tex
	new file:   tex/d2q5/d2q5m5_matB.tex
	new file:   tex/d2q5/d2q5m5_matB_per_els.tex
	new file:   tex/d2q5/d2q5m5_matK.tex
	new file:   tex/d2q5/d2q5m5_matK_per_els.tex
	new file:   tex/d2q5/d2q5m6_feq.tex
	new file:   tex/d2q5/d2q5m6_matA.tex
	new file:   tex/d2q5/d2q5m6_matA_per_els.tex
	new file:   tex/d2q5/d2q5m6_matB.tex
	new file:   tex/d2q5/d2q5m6_matB_per_els.tex
	new file:   tex/d2q5/d2q5m6_matK.tex
	new file:   tex/d2q5/d2q5m6_matK_per_els.tex
	new file:   tex/d2q5/supp_d2q5_defs.tex
	new file:   tex/d2q5/supp_d2q5m0_defs.tex
	new file:   tex/d2q5/supp_d2q5m1_defs.tex
	new file:   tex/d2q5/supp_d2q5m2_defs.tex
	new file:   tex/d2q5/supp_d2q5m3_defs.tex
	new file:   tex/d2q5/supp_d2q5m4_defs.tex
	new file:   tex/d2q5/supp_d2q5m5_defs.tex
	new file:   tex/d2q5/supp_d2q5m6_defs.tex
	new file:   tex/d3q7/d3q7_C.tex
	new file:   tex/d3q7/d3q7_defs.tex
	new file:   tex/d3q7/d3q7_matM.tex
	new file:   tex/d3q7/d3q7_matM_per_els.tex
	new file:   tex/d3q7/d3q7m0_feq.tex
	new file:   tex/d3q7/d3q7m0_matA.tex
	new file:   tex/d3q7/d3q7m0_matA_per_els.tex
	new file:   tex/d3q7/d3q7m0_matB.tex
	new file:   tex/d3q7/d3q7m0_matB_per_els.tex
	new file:   tex/d3q7/d3q7m0_matM.tex
	new file:   tex/d3q7/d3q7m0_matM_per_els.tex
	new file:   tex/d3q7/d3q7m1_feq.tex
	new file:   tex/d3q7/d3q7m1_matA.tex
	new file:   tex/d3q7/d3q7m1_matA_per_els.tex
	new file:   tex/d3q7/d3q7m1_matB.tex
	new file:   tex/d3q7/d3q7m1_matB_per_els.tex
	new file:   tex/d3q7/d3q7m1_matM.tex
	new file:   tex/d3q7/d3q7m1_matM_per_els.tex
	new file:   tex/d3q7/d3q7m2_feq.tex
	new file:   tex/d3q7/d3q7m2_matA.tex
	new file:   tex/d3q7/d3q7m2_matA_per_els.tex
	new file:   tex/d3q7/d3q7m2_matB.tex
	new file:   tex/d3q7/d3q7m2_matB_per_els.tex
	new file:   tex/d3q7/d3q7m2_matM.tex
	new file:   tex/d3q7/d3q7m2_matM_per_els.tex
	new file:   tex/d3q7/d3q7m3_feq.tex
	new file:   tex/d3q7/d3q7m3_matA.tex
	new file:   tex/d3q7/d3q7m3_matA_per_els.tex
	new file:   tex/d3q7/d3q7m3_matB.tex
	new file:   tex/d3q7/d3q7m3_matB_per_els.tex
	new file:   tex/d3q7/d3q7m3_matM.tex
	new file:   tex/d3q7/d3q7m3_matM_per_els.tex
	new file:   tex/d3q7/d3q7m4_feq.tex
	new file:   tex/d3q7/d3q7m4_matA.tex
	new file:   tex/d3q7/d3q7m4_matA_per_els.tex
	new file:   tex/d3q7/d3q7m4_matB.tex
	new file:   tex/d3q7/d3q7m4_matB_per_els.tex
	new file:   tex/d3q7/d3q7m4_matM.tex
	new file:   tex/d3q7/d3q7m4_matM_per_els.tex
	new file:   tex/d3q7/d3q7m5_feq.tex
	new file:   tex/d3q7/d3q7m5_matA.tex
	new file:   tex/d3q7/d3q7m5_matA_per_els.tex
	new file:   tex/d3q7/d3q7m5_matB.tex
	new file:   tex/d3q7/d3q7m5_matB_per_els.tex
	new file:   tex/d3q7/d3q7m5_matK.tex
	new file:   tex/d3q7/d3q7m5_matK_per_els.tex
	new file:   tex/d3q7/d3q7m5_matM.tex
	new file:   tex/d3q7/d3q7m5_matM_per_els.tex
	new file:   tex/d3q7/d3q7m6_feq.tex
	new file:   tex/d3q7/d3q7m6_matA.tex
	new file:   tex/d3q7/d3q7m6_matA_per_els.tex
	new file:   tex/d3q7/d3q7m6_matB.tex
	new file:   tex/d3q7/d3q7m6_matB_per_els.tex
	new file:   tex/d3q7/d3q7m6_matK.tex
	new file:   tex/d3q7/d3q7m6_matK_per_els.tex
	new file:   tex/d3q7/d3q7m6_matM.tex
	new file:   tex/d3q7/d3q7m6_matM_per_els.tex
	new file:   tex/d3q7/supp_d3q7_defs.tex
	new file:   tex/d3q7/supp_d3q7m0_defs.tex
	new file:   tex/d3q7/supp_d3q7m1_defs.tex
	new file:   tex/d3q7/supp_d3q7m2_defs.tex
	new file:   tex/d3q7/supp_d3q7m3_defs.tex
	new file:   tex/d3q7/supp_d3q7m4_defs.tex
	new file:   tex/d3q7/supp_d3q7m5_defs.tex
	new file:   tex/d3q7/supp_d3q7m6_defs.tex
	new file:   tex/header.tex
	new file:   tex/supp_d1q3c_efde.pdf
	new file:   tex/supp_d1q3c_efde.tex
	new file:   tex/supp_d1q3c_epde.pdf
	new file:   tex/supp_d1q3c_epde.tex
	new file:   tex/supp_d1q3n_efde.tex
	new file:   tex/supp_d1q3n_epde.tex
	new file:   tex/supp_d2q5c_efde.pdf
	new file:   tex/supp_d2q5c_efde.tex
	new file:   tex/supp_d2q5c_epde.pdf
	new file:   tex/supp_d2q5c_epde.tex
	new file:   tex/supp_d2q5n_efde.tex
	new file:   tex/supp_d2q5n_epde.tex
	new file:   tex/supp_d3q7c_efde.tex
	new file:   tex/supp_d3q7c_epde.tex
	new file:   tex/supp_d3q7n_efde.tex
	new file:   tex/supp_d3q7n_epde.tex
parent 9bbc4789

coefs.h

+2 −0
Original line number Diff line number Diff line
@@ -25,6 +25,8 @@ struct COEFS 
	std::vector<SC> v;

	std::vector<ex> basis;

	std::vector<ex> basis_tex;

	std::vector<int> pow_dl;

	std::vector<int> pow_dt;

	int order=6; // default

	void fill(ex &e);

	bool fill_recursive(ex &e, SC &sc, int level);

core.cpp

+40 −4
Original line number Diff line number Diff line
@@ -153,6 +153,8 @@ void CORE::fill_sym_vels_to_basis(COEFS &CF) 
	{

		CF.basis.push_back(sym_vel[d]);

		CF.basis_tex.push_back(sym_vel[d]);

		CF.pow_dl.push_back(0);

		CF.pow_dt.push_back(0);

	}



	if (use_tay_vels)
@@ -165,6 +167,9 @@ void CORE::fill_sym_vels_to_basis(COEFS &CF) 
	{

		for (int d=0;d<DIM;d++)

		{

			CF.pow_dl.push_back(tay_vels[d]->get_pow_dl(k-j,j1,j2,j3));

			CF.pow_dt.push_back(tay_vels[d]->get_pow_dt(k-j,j1,j2,j3));

			

			CF.basis.push_back(tay_vels[d]->get_tay(k-j,j1,j2,j3));

			CF.basis_tex.push_back(tay_vels[d]->get_tay_tex(k-j,j1,j2,j3));

		}
@@ -251,7 +256,7 @@ void CORE::allocate() 
	{

		tay_vels = new TAYLOR*[DIM];

		if (verbosity==1) printf("allocate tay_vels\n");

		for (int d=0;d<DIM;d++) tay_vels[d] = new TAYLOR(DIM,Q,ORDER,_index,Nindex);

		for (int d=0;d<DIM;d++) tay_vels[d] = new TAYLOR(DIM,Q,ORDER,_index,Nindex,sym_dl,sym_dt);

		if (DIM>=1) tay_vels[0]->allocate("u", sym_vel[0]);

		if (DIM>=2) tay_vels[1]->allocate("v", sym_vel[1]);

		if (DIM>=3) tay_vels[2]->allocate("w", sym_vel[2]);
@@ -350,7 +355,7 @@ void CORE::allocate() 
	{

		if (use_tay_ms)

		{

			tay_ms[q] = new TAYLOR(DIM,Q,ORDER,_index,Nindex);

			tay_ms[q] = new TAYLOR(DIM,Q,ORDER,_index,Nindex,sym_dl,sym_dt);

			lbm.moment_symbol(desc2,q,Q);

			lbm.moment_symbol_number(desc3,q,Q);

			if (strlen(desc3)>0)
@@ -361,7 +366,7 @@ void CORE::allocate() 
		}

		if (use_tay_meqs)

		{

			tay_meqs[q] = new TAYLOR(DIM,Q,ORDER,_index,Nindex);

			tay_meqs[q] = new TAYLOR(DIM,Q,ORDER,_index,Nindex,sym_dl,sym_dt);

			sprintf(comp,"{{\\tilde m}^{eq}_{%d}}",q); 

			tay_meqs[q]->allocate(comp,0);

		}
@@ -1079,8 +1084,25 @@ void CORE::export_EPDE_coefs_bigprint(ex &expr, ofstream &fout, char* str) 
			}

			if (none) fout << "1";

			fout << "\\right]} }";

			

						// put dl^pow dt^pow ... count all powers

			int pow_dl=0;

			int pow_dt=0;



			if (order>0) 

			{

				for (int b=0;b<CF.basis.size();b++)

				if (CF.v[i].exp[b] > 0)

				{

					pow_dl+=CF.pow_dl[b]*CF.v[i].exp[b];

					pow_dt+=CF.pow_dt[b]*CF.v[i].exp[b];

				}

				// put dl^pow dt^pow

				if (pow_dt==1) fout << "\\dt ";

				if (pow_dt>1)  fout << "\\dt^{"<< pow_dt <<"} ";

				if (pow_dl==1) fout << "\\dl ";

				if (pow_dl>1)  fout << "\\dl^{"<< pow_dl <<"} ";



				for (int b=0;b<CF.basis.size();b++)

				{

					if (CF.v[i].exp[b] != 0)
@@ -1134,6 +1156,12 @@ void CORE::export_EPDE_coefs_bigprint(ex &expr, ofstream &fout, char* str) 
			if (CF2.v[i2].order()==or2)

				outex += CF2.v[i2].coef * CF2.v[i2].base;



			outex = expand(outex);



			// eliminate \dl and \dt from the coefficient

			outex = subs(outex, sym_dl==1);

			outex = subs(outex, sym_dt==1);



			ostringstream s1;

			s1 << latex;

			if (factor_coefs)
@@ -1196,6 +1224,14 @@ void CORE::export_EPDE_coefs(int pde) 
			if (k-j>0) { if (k-j==1) sprintf(str_dt,"\\dt"); else sprintf(str_dt,"\\dt^{%d}",k-j); } else sprintf(str_dt,"");

			if (j>0) { if (j==1) sprintf(str_dl,"\\dl"); else sprintf(str_dl,"\\dl^{%d}",j); } else sprintf(str_dl,"");

			

			// put dl^pow dt^pow

			int pow_dl=tay_ms[q]->get_pow_dl(k-j,j1,j2,j3);

			int pow_dt=tay_ms[q]->get_pow_dt(k-j,j1,j2,j3);

			if (pow_dt==1) gout << "\\dt ";

			if (pow_dt>1) gout << "\\dt^{"<< pow_dt <<"} ";

			if (pow_dl==1) gout << "\\dl ";

			if (pow_dl>1) gout << "\\dl^{"<< pow_dl <<"} ";

			

			gout << tay_ms[q]->get_tay(k-j,j1,j2,j3) << endl;

			last_c_nonzero=count;

			count++;
@@ -1216,7 +1252,7 @@ void CORE::export_EPDE_coefs(int pde) 
	for (int j1=0;j1<=k;j1++)

	for (int q=0;q<Q;q++)

	if (j1+j2+j3==j)

		if (coefs_epde[pde][k-j][indO(j1,j2,j3)]!=0)

		if (coefs_epde[pde][q][k-j][indO(j1,j2,j3)]!=0)

			max++;

	

	count=0;

core.h

+3 −0
Original line number Diff line number Diff line
@@ -66,6 +66,9 @@ struct CORE : DEF 
	// central values:

	ex *sym_m, *sym_meq, *sym_f, *sym_feq;



	ex sym_dl=symbol("\\dl");

	ex sym_dt=symbol("\\dt");



	TAYLOR **tay_vels, **tay_ms, **tay_meqs;



	matrix const_feq, const_sub_feq, const_meq;

lbm.cpp

+533 −566

File changed.

Preview size limit exceeded, changes collapsed.

main.cpp

+56 −15
Original line number Diff line number Diff line
@@ -29,22 +29,20 @@ 
#include "defs.h"

#include "core.h"



int main(int argc, char* argv[]) 

void execute(

	int DIM = 3,			// dimension

	int Q = 7,			// number of discrete velocities

	int model = SRT,		// collision model

	bool constvel = true,		// velocity field is constant in time and space (for ADE)

	bool bigprint = false,		// use a more detailed decomposition of exported coefficients

	bool verbosity = true,

	bool EFDEprint = true,		// export EFDE

	bool EPDEprint = true,		// export EPDE

	bool factor_coefs = false,	// use GiNaC::factor() when exporting coefficients to latex

	int feqorder = 2, 		// -1 = symbolic, 1 = linear, 2 = quadratic

	int ORDER = 2			// Taylor expansion order

)

{

	int DIM = 3;			// dimension

	int Q = 7;			// number of discrete velocities

	int model = SRT;		// collision model

	int ORDER = 2;			// Taylor expansion order



	bool verbosity = true;

	bool bigprint = true;		// use a more detailed decomposition of exported coefficients

	bool constvel = true;		// velocity field is constant in time and space (for ADE)

	bool EFDEprint = true;		// export EFDE

	bool EPDEprint = true;		// export EPDE

	bool factor_coefs = false;	// use GiNaC::factor() when exporting coefficients to latex

	int feqorder = -1; //-1 = symbolic





	printf("d=%d q=%d model=%d feqorder=%d constvel=%s bigprint=%d\n", DIM, Q, model, feqorder, (constvel)?"yes":"no", bigprint);

	fflush(stdout);
@@ -89,6 +87,49 @@ int main(int argc, char* argv[]) 
	if (verbosity) { printf("EPDE coefs computed\n"); fflush(stdout); }



	if (EPDEprint) for (int d=0;d<Q;d++) core.export_EPDE_coefs(d);

}



void d1q3(bool constvel=true)

{

	execute(1,3,SRT,constvel);

	execute(1,3,MRT1,constvel);

	execute(1,3,CLBM1,constvel);

}



void d2q5(bool constvel=true)

{

	bool bigprint=!constvel;

	execute(2,5,SRT,constvel,bigprint);

	execute(2,5,MRT1,constvel,bigprint);

	execute(2,5,MRT2,constvel,bigprint);

	execute(2,5,MRT1_OG,constvel,bigprint);

	execute(2,5,MRT2_OG,constvel,bigprint);

	execute(2,5,CLBM1,constvel,bigprint);

	execute(2,5,CLBM2,constvel,bigprint);

}



void d3q7(bool constvel=true)

{

	bool bigprint=!constvel;

	execute(3,7,SRT,constvel,bigprint);

	execute(3,7,MRT1,constvel,bigprint);

	execute(3,7,MRT2,constvel,bigprint);

	execute(3,7,MRT1_OG,constvel,bigprint);

	execute(3,7,MRT2_OG,constvel,bigprint);

	execute(3,7,CLBM1,constvel,bigprint);

	execute(3,7,CLBM2,constvel,bigprint);

}



int main(int argc, char* argv[]) 

{

	d1q3(true);

//	d1q3(false);



	d2q5(true);

//	d2q5(false);



//	d3q7(true);

//	d3q7(false);

	return 0;

}