Simulação do Modelo de Gray-Scott - Histórico de revisão

Wallec em 21h54min de 4 de março de 2022

2022-03-04T21:54:12Z

← Edição anterior		Edição das 18h54min de 4 de março de 2022
Linha 196:		Linha 196:
	printf("set palette gray\n");		printf("set palette gray\n");
	printf("set pm3d map\n");		printf("set pm3d map\n");
			printf("set cbrange [0:1]\n");
	printf("splot \"-\" matrix w image\n");		printf("splot \"-\" matrix w image\n");

Wallec em 19h19min de 4 de março de 2022

2022-03-04T19:19:02Z

← Edição anterior		Edição das 16h19min de 4 de março de 2022
Linha 201:		Linha 201:
	{		{
	for(int j=0; j<L; j++)		for(int j=0; j<L; j++)
	printf("%lf ", (arr+iL+j));		printf("%lf ", (arr+i+jL));

	printf("\n");		printf("\n");

Wallec em 19h16min de 4 de março de 2022

2022-03-04T19:16:04Z

@@ Linha 127: / Linha 127: @@
 /*Atualização das matrizes. Esse é o núcleo do programa. Integram-se as
-equações do modelo com FTCS. São usadas as seguintes condições de fronteira:
+equações do modelo com FTCS. São usadas as seguintes condições de fronteira
 periódicas:
 * Se i=0 (primeira linha da matriz), o vizinho da esquerda está na última
@@ Linha 154: / Linha 141: @@
 * Se i=L-1 (última linha da matriz, o vizinho da direita está na primeira
 linha da matriz, com mesmo j).
 */
@@ Linha 170: / Linha 153: @@
      for(j=0; j<L; j++)
+     {
-       /*Condições totalmente periódicas:
+       /*Condições de fronteira periódicas:*/
        uC = *(u+i*L+j);
@@ Linha 183: / Linha 166: @@
        vL = (i==0)? *(v+L2-L+j) : *(v+(i-1)*L+j);
        vD = (j==0)? *(v+i*L+(L-1)) : *(v+i*L+(j-1));
        /*Laplacianos discretizados*/

Wallec em 18h29min de 2 de março de 2022

2022-03-02T18:29:42Z

← Edição anterior		Edição das 15h29min de 2 de março de 2022
Linha 127:		Linha 127:

	/*Atualização das matrizes. Esse é o núcleo do programa. Integram-se as		/*Atualização das matrizes. Esse é o núcleo do programa. Integram-se as
	equações do modelo com FTCS. São usadas condições de fronteira ~~conforme o Sayama~~:		equações do modelo com FTCS. São usadas as seguintes condições de fronteira:

	* Se i=0 (primeira linha da matriz), o vizinho da esquerda é o vizinho da direita;		* Se i=0 (primeira linha da matriz), o vizinho da esquerda é o vizinho da direita;

Wallec em 18h28min de 2 de março de 2022

2022-03-02T18:28:05Z

← Edição anterior		Edição das 15h28min de 2 de março de 2022
Linha 175:		Linha 175:
	uU = (u+iL+(j+1)%L);		uU = (u+iL+(j+1)%L);
	uR = (u+((i+1)%L)L+j);		uR = (u+((i+1)%L)L+j);
	uL = (i==0)? (u+L2-L+j) : (u+(i-1)L+j);		uL = (i==0)? (u+L2-L+j) : (u+(i-1)*L+j);
	uD = (j==0)? (u+iL+(L-1)) : (u+iL+(j-1));		uD = (j==0)? (u+iL+(L-1)) : (u+iL+(j-1));

Linha 186:		Linha 186:
	*/		*/

	/Condições ~~do Sayama~~/		/Condições usadas/
	uC = (u+iL+j);		uC = (u+iL+j);
	uU = (u+iL+(j+1)%L);		uU = (u+iL+(j+1)%L);

Wallec em 20h37min de 26 de fevereiro de 2022

2022-02-26T20:37:34Z

← Edição anterior		Edição das 17h37min de 26 de fevereiro de 2022
Linha 35:		Linha 35:
	const int L2=4761;		const int L2=4761;
	const double dh=0.014;		const double dh=0.014;
	const double dt=0.01;		const double dt=0.02;

	double tmax;		double tmax;
Linha 64:		Linha 64:
	{		{
	update(u, v, n_u, n_v);		update(u, v, n_u, n_v);

	~~/Se usado -DGNU na compilação, executa o código abaixo/~~
	#ifdef GNU		#ifdef GNU
	if(count % 50 == 0)		if(count % 50 == 0)
Linha 78:		Linha 77:
	}		}

	~~/Se usado -DGNU na compilação, executa o código abaixo/~~
	#ifdef GNU		#ifdef GNU
	gnuplot(u);		gnuplot(u);
	printf("pause 5\n");		printf("pause 5\n");
	#endif		#endif

	~~/Gera a saída final (por enquanto, apenas de u), no tempo ~tmax, no arquivo./~~

	sprintf(output,"gs_f_%.3f_k_%.3f_tmax_%.0f.dat",f, k, tmax);		sprintf(output,"gs_f_%.3f_k_%.3f_tmax_%.0f.dat",f, k, tmax);
Linha 120:		Linha 116:
	{		{
	(u+iL+j) = (i*L+j == (L2-1)/2)? 0.0 : 1.0;		(u+iL+j) = (i*L+j == (L2-1)/2)? 0.0 : 1.0;
	(v+iL+j) = (i*L+j ~~== (L2-1)/2~~)~~? 1.0 : 0.0~~;		(v+iL+j) = 1.0-(u+iL+j);
	}		}

Linha 131:		Linha 127:

	/*Atualização das matrizes. Esse é o núcleo do programa. Integram-se as		/*Atualização das matrizes. Esse é o núcleo do programa. Integram-se as
	equações do modelo com FTCS. São usadas condições de ~~contorno~~ periódicas.		equações do modelo com FTCS. São usadas condições de fronteira conforme o Sayama:

			* Se i=0 (primeira linha da matriz), o vizinho da esquerda é o vizinho da direita;

			* Se j=0 (primeira coluna da matriz), o vizinho de baixo é o vizinho de cima;

			* Se j=L-1 (última coluna da matriz), o vizinho de cima está na primeira
			coluna da matriz, com mesmo i;

			* Se i=L-1 (última linha da matriz, o vizinho da direita está na primeira
			linha da matriz, com mesmo j).

			*********************************************************************************

			Também testamos várias simulações com condições completamente periódicas:

	* Se i=0 (primeira linha da matriz), o vizinho da esquerda está na última		* Se i=0 (primeira linha da matriz), o vizinho da esquerda está na última
Linha 144:		Linha 154:
	* Se i=L-1 (última linha da matriz, o vizinho da direita está na primeira		* Se i=L-1 (última linha da matriz, o vizinho da direita está na primeira
	linha da matriz, com mesmo j).		linha da matriz, com mesmo j).

			Os comportamentos gerais da formação de padrões são os mesmos, mas as figuras
			diferem, principalmente para tempos longos, em que a difusão atinge as
			fronteiras do grid.
	*/		*/

Linha 156:		Linha 170:
	for(j=0; j<L; j++)		for(j=0; j<L; j++)
	{		{
			/*Condições totalmente periódicas:

	uC = (u+iL+j);		uC = (u+iL+j);
	uR = ~~(i==L-1)?~~ (u+~~j) :~~ (u+(i+1)*L+j);		uU = (u+iL+(j+1)%L);
	uU = ~~(j==L-1)?~~ (u+iL) ~~: (u+i~~L+j+1);		uR = (u+((i+1)%L)L+j);
	uL = (i==0)? (u+L2-L+j) : (u+(i-1)*L+j);		uL = (i==0)? (u+L2-L+j) : (u+(i-1)L+j);
	uD = (j==0)? (u+iL+(L-1)) : (u+iL+(j-1));		uD = (j==0)? (u+iL+(L-1)) : (u+iL+(j-1));

	vC = (v+iL+j);		vC = (v+iL+j);
	~~vR = (i~~=~~=L-1)?~~ (v+~~j) :~~ (v+(i+1)*L+j);		vU = (v+iL+(j+1)%L);
	~~vU = (j~~=~~=L-1)?~~ (v+iL) ~~: (v+i~~L+j+1);		vR = (v+((i+1)%L)L+j);
	vL = (i==0)? (v+L2-L+j) : (v+(i-1)*L+j);		vL = (i==0)? (v+L2-L+j) : (v+(i-1)*L+j);
	vD = (j==0)? (v+iL+(L-1)) : (v+iL+(j-1));		vD = (j==0)? (v+iL+(L-1)) : (v+iL+(j-1));

			*/

			/Condições do Sayama/
			uC = (u+iL+j);
			uU = (u+iL+(j+1)%L);
			uR = (u+((i+1)%L)L+j);
			uL = (u+abs((i-1)%L)L+j);
			uD = (u+iL+abs((j-1)%L));

			vC = (v+iL+j);
			vU = (v+iL+(j+1)%L);
			vR = (v+((i+1)%L)L+j);
			vL = (v+abs((i-1)%L)L+j);
			vD = (v+iL+abs((j-1)%L));

	/Laplacianos discretizados/		/Laplacianos discretizados/
	uLap = (uR + uL + uU + uD - 4uC)/(dhdh);		uLap = (uR + uL + uU + uD - 4uC)/(dhdh);
Linha 185:		Linha 216:
	}		}

	/Visualização no gnuplot. Autoexplicativo./		/*Visualização no gnuplot. Autoexplicativo. Este módulo é praticamente
			igual ao criado pelo grupo do modelo de Turing. Apenas foi acrescentado
			o comando set palette gray*/
	void gnuplot(double *arr)		void gnuplot(double *arr)
	{		{
Linha 193:		Linha 226:
	printf("set size square\n");		printf("set size square\n");
	printf("set key off\n");		printf("set key off\n");
	printf("set palette ~~grey~~\n");		printf("set palette gray\n");
	printf("~~plot~~ \"-\" matrix w image\n");		printf("set pm3d map\n");
			printf("splot \"-\" matrix w image\n");

	for(int i=0; i<L; i++)		for(int i=0; i<L; i++)

Wallec em 01h24min de 25 de fevereiro de 2022

2022-02-25T01:24:35Z

← Edição anterior		Edição das 22h24min de 24 de fevereiro de 2022
Linha 64:		Linha 64:
	{		{
	update(u, v, n_u, n_v);		update(u, v, n_u, n_v);

			/Se usado -DGNU na compilação, executa o código abaixo/
	#ifdef GNU		#ifdef GNU
	if(count % 50 == 0)		if(count % 50 == 0)
Linha 77:		Linha 78:
	}		}

			/Se usado -DGNU na compilação, executa o código abaixo/
	#ifdef GNU		#ifdef GNU
	gnuplot(u);		gnuplot(u);
	printf("pause 5\n");		printf("pause 5\n");
	#endif		#endif

			/Gera a saída final (por enquanto, apenas de u), no tempo ~tmax, no arquivo./

	sprintf(output,"gs_f_%.3f_k_%.3f_tmax_%.0f.dat",f, k, tmax);		sprintf(output,"gs_f_%.3f_k_%.3f_tmax_%.0f.dat",f, k, tmax);
Linha 164:		Linha 168:
	vD = (j==0)? (v+iL+(L-1)) : (v+iL+(j-1));		vD = (j==0)? (v+iL+(L-1)) : (v+iL+(j-1));

	//Laplacianos discretizados		/Laplacianos discretizados/
	uLap = (uR + uL + uU + uD - 4uC)/(dhdh);		uLap = (uR + uL + uU + uD - 4uC)/(dhdh);
	vLap = (vR + vL + vU + vD - 4vC)/(dhdh);		vLap = (vR + vL + vU + vD - 4vC)/(dhdh);

	//Equações do modelo discretizadas		/Equações do modelo discretizadas/
	(n_u+iL+j) = uC + (-uCvCvC + f(1-uC) + DuuLap)*dt;		(n_u+iL+j) = uC + (-uCvCvC + f(1-uC) + DuuLap)*dt;
	(n_v+iL+j) = vC + (uCvCvC - (f+k)vC + DvvLap)*dt;		(n_v+iL+j) = vC + (uCvCvC - (f+k)vC + DvvLap)*dt;

Wallec em 01h15min de 25 de fevereiro de 2022

2022-02-25T01:15:46Z

← Edição anterior		Edição das 22h15min de 24 de fevereiro de 2022
Linha 204:		Linha 204:
	}		}

	<source/>		</source>

Wallec em 01h12min de 25 de fevereiro de 2022

2022-02-25T01:12:14Z

← Edição anterior		Edição das 22h12min de 24 de fevereiro de 2022
Linha 1:		Linha 1:
	<lang = "c">		<source lang = "c">

	/********************************************************		/********************************************************
Linha 203:		Linha 203:

	}		}

			<source/>

Wallec: Criou página com ' /******************************************************** * Modelo de Gray-Scott * *************************************************/ /...'

2022-02-25T01:10:31Z

Criou página com '<lang = "c"> ******************************************************* * Modelo de Gray-Scott * ******************************************************: /******...'

Página nova

<lang = "c">

/********************************************************
* Modelo de Gray-Scott *
*******************************************************/
/********************************************************
* Grid de 69 x 69 (L x L) *
********************************************************
* Parâmetros (na ordem): *
* tmax F k Du Dv *
* *
* Exemplo: *
* ./a.out 1000 0.015 0.055 0.00002 0.00001 *
* ./a.out tmax F k Du Dv *
* *
* Saída para gnuplot: *
********************************************************
* compilar com: *
* gcc -DGNU gray-scott.c -lm *
********************************************************
* Executar com pipe: *
* ./a.out PARAMETROS | gnuplot *
* *
********************************************************
* O programa sempre gera um arquivo com o *
* resultado da simulação para tmax. Exemplo: *
* gs_f_0.015_k_0.055_tmax_1000.dat *
********************************************************/

#include<stdio.h>
#include<stdlib.h>
#include<math.h>

const int L=69;
const int L2=4761;
const double dh=0.014;
const double dt=0.01;

double tmax;
double Du, Dv;
double f, k;

void init(double *, double *, char *[]);
void update(double *, double *, double *, double *);
void gnuplot(double *);

int main(int argc, char *argv[ ])
{
FILE *fp;
char output[100];
double t = 0;
double *u, *v;
double *n_u, *n_v;
int count = 0;

u = (double*)calloc(L2, sizeof(double));
v = (double*)calloc(L2, sizeof(double));
n_u = (double*)calloc(L2, sizeof(double));
n_v = (double*)calloc(L2, sizeof(double));

init(u, v, argv);

while(t < tmax)
{
update(u, v, n_u, n_v);

#ifdef GNU
if(count % 50 == 0)
{
printf("set title 'Tempo %.2lf'\n", t);
gnuplot(u);
}
#endif

t += dt;
count += 1;
}

#ifdef GNU
gnuplot(u);
printf("pause 5\n");
#endif

sprintf(output,"gs_f_%.3f_k_%.3f_tmax_%.0f.dat",f, k, tmax);
fp = fopen(output, "w");

for(int i=0; i<L; i++)
{
for(int j=0; j<L; j++)
fprintf(fp, "%lf ", *(u+i*L+j));

fprintf(fp, "\n");
}

fclose(fp);

free(u);
free(v);
free(n_u);
free(n_v);

return 0;
}

/*INICIALIZAÇÃO:
Escolhendo o ponto do meio do grid para tirar o sistema do equilíbrio,
os demais pontos estão no equilíbrio trivial (u, v) = (1, 0).
Pressupõe-se que L seja ímpar. Se L for par, o ponto inicialmente fora
do equilíbrio não será central.
*/

void init(double *u, double *v, char *argv[ ])
{
for(int i=0; i<L; i++)
for(int j=0; j<L; j++)
{
*(u+i*L+j) = (i*L+j == (L2-1)/2)? 0.0 : 1.0;
*(v+i*L+j) = (i*L+j == (L2-1)/2)? 1.0 : 0.0;
}

tmax = atof(argv[1]);
f = atof(argv[2]);
k = atof(argv[3]);
Du = atof(argv[4]);
Dv = atof(argv[5]);
}

/*Atualização das matrizes. Esse é o núcleo do programa. Integram-se as
equações do modelo com FTCS. São usadas condições de contorno periódicas.

* Se i=0 (primeira linha da matriz), o vizinho da esquerda está na última
linha da matriz, com mesmo j;

* Se j=0 (primeira coluna da matriz), o vizinho de baixo está na última
coluna da matriz, com mesmo i;

* Se j=L-1 (última coluna da matriz), o vizinho de cima está na primeira
coluna da matriz, com mesmo i

* Se i=L-1 (última linha da matriz, o vizinho da direita está na primeira
linha da matriz, com mesmo j).
*/

void update(double *u, double *v, double *n_u, double *n_v)
{
int i, j;
double uLap, vLap;
double uC, uR, uL, uU, uD;
double vC, vR, vL, vU, vD;

for(i=0; i<L; i++)
for(j=0; j<L; j++)
{
uC = *(u+i*L+j);
uR = (i==L-1)? *(u+j) : *(u+(i+1)*L+j);
uU = (j==L-1)? *(u+i*L) : *(u+i*L+j+1);
uL = (i==0)? *(u+L2-L+j) : *(u+(i-1)*L+j);
uD = (j==0)? *(u+i*L+(L-1)) : *(u+i*L+(j-1));

vC = *(v+i*L+j);
vR = (i==L-1)? *(v+j) : *(v+(i+1)*L+j);
vU = (j==L-1)? *(v+i*L) : *(v+i*L+j+1);
vL = (i==0)? *(v+L2-L+j) : *(v+(i-1)*L+j);
vD = (j==0)? *(v+i*L+(L-1)) : *(v+i*L+(j-1));

//Laplacianos discretizados
uLap = (uR + uL + uU + uD - 4*uC)/(dh*dh);
vLap = (vR + vL + vU + vD - 4*vC)/(dh*dh);

//Equações do modelo discretizadas
*(n_u+i*L+j) = uC + (-uC*vC*vC + f*(1-uC) + Du*uLap)*dt;
*(n_v+i*L+j) = vC + (uC*vC*vC - (f+k)*vC + Dv*vLap)*dt;
}

for(i=0; i<L; i++)
for(j=0; j<L; j++)
{
*(u+i*L+j) = *(n_u+i*L+j);
*(v+i*L+j) = *(n_v+i*L+j);
}
}

/*Visualização no gnuplot. Autoexplicativo.*/
void gnuplot(double *arr)
{
printf("set xrange [0:%d]\nset yrange [0:%d]\n", (L-1), (L-1));
printf("unset xtics\nunset ytics\n");
printf("unset colorbox\n");
printf("set size square\n");
printf("set key off\n");
printf("set palette grey\n");
printf("plot \"-\" matrix w image\n");

for(int i=0; i<L; i++)
{
for(int j=0; j<L; j++)
printf("%lf ", *(arr+i*L+j));

printf("\n");
}

printf("\n\n\ne\n\n");

}