Condução de Calor em 2D - Código somente funciona para determinados tamanhos de matriz

KingOfCoke
(usa Outra)

Enviado em 28/10/2016 - 18:06h

Olá, boa tarde.

Sou aluno de Eng Civil (só tive o básico de programação em C) e estou fazendo um trabalho de faculdade sobre condução de calor em 2D utilizando diferenças finitas e SOR para a matéria de Cálculo Numérico.

Aparentemente o código funciona, porém somente para alguns tamanhos de matrizes. Defini NX e NY como constantes, já que vou ter de fazer o cálculo para diferentes tamanhos malhas. Eis o problema: em alguns tamanhos funcionam (testei NX = NY = 3, 4, 6, 8, 16, 18), mas em outros o erro da igual à -1.#IND000000 e/ou os valores T[i][j] estouram para um valor absurdo. Nesse código em específico, estou precisando para NX = NY = 10.

Alguém poderia me ajudar? Ou me dar um norte do por quê isso ocorre?

Lembrando que em outros tamanhos de função o código funciona. Testei 18x18 e comparando com a solução exata que é dada, os valores são aceitáveis.

Obrigado!

#include <stdio.h>

#include <stdlib.h>

#include <conio.h>

#include <math.h>

#define NX 10

#define NY 10

#define K 2

#define HX 0.111111111

#define HY 0.111111111





typedef struct

{

	float w;

	int n;

	double erro;

	double T[NX][NY];

} dataFromSOR;





typedef struct

{

	double F[NX][NY];

} FMatrix;





void interface () 

{

	printf("\n  Validacao de um modelo conhecido de conducao de calor em 2D . \n  Sendo a equacao de T dada e f(x,y) encontrado atraves de tal .\n\n  Apresentacao de solucao em arquivo .txt para w otimo .\n\n  Por favor, aguarde...\n\n");

}





FMatrix findF ()

{

	double xi, yj;

	int i, j;

	FMatrix FM;

	

	for (j=0 ; j<NY; j++)

	{

		for (i=0 ; i<NX; i++)

		{

			xi = i*HX;

			yj = j*HY;

		

			FM.F[i][j] = 200*(-K)*( xi*xi + yj*yj - xi - yj);

	

		}

	}

	return FM;

}



dataFromSOR SOR (float w, FMatrix FM, double a, double b, double c, double d, double e, double f)

{

	double T[NX][NY], D[NX][NY], maxdif, maxvet, erro;

	int i, n, j;

	dataFromSOR data;



		n = 0;

		erro = 300; 

	

		

		j=i=0;

		while (j<NY+1)

		{

		 	while (i<NX+1) 

			 {

			 	

		 		D[i][j] = 0;               

		 		T[i][j] = 0;

				 i++;           

			 }

			 j++;

		} 

			

		

		while (erro > 0.0000001)

		{

				

			maxdif = 0;

			maxvet = 0;

			

			for (j=2; j<NY; j++)

			{

				

				for (i=2; i<NX; i++)

				{

										

					D[i][j] = T[i][j];

									

					T[i][j] = T[i][j]*(1 - w) + (w/e)*( f*FM.F[i][j] - (a*T[i+1][j] + b*T[i-1][j] + c*T[i][j+1] + d*T[i][j-1]));

					D[i][j] = fabs(D[i][j] - T[i][j]);

					data.T[i][j] = T[i][j];

					

					printf("F[%d][%d] = %2.1f \t", i, j, T[i][j]);

					

					if (fabs(T[i][j]) > maxvet)

					{

						maxvet = T[i][j];

					}

					

//					printf("maxvet[%d][%d] = %2.1f \t", i, j, maxvet);

					

					if(D[i][j] > maxdif)

					{

						maxdif = D[i][j];

					}

//					printf("maxdif[%d][%d] = %2.1f \t", i, j, maxdif);

				}

				printf("\n");

			}

			

		

			erro = maxdif/maxvet;        

			n++;                             

		}

		

		

		data.w = w;

		data.n = n;

		data.erro = erro;

		

		printf("w = %2.2f     n = %d     Erro = %2.10f \n" , w, n, erro);

		

		return data;

}





void writeW (dataFromSOR data)

{





	FILE *dados;

	dados = fopen ("dados8x8.txt","a");

		

	if (dados == NULL)

	{

		printf ("Houve um erro ao criar o arquivo dados.\n");

	}

	

	fprintf(dados, "%2.2f \t %d \t %2.10f \n", data.w, data.n, data.erro);

				

	if(fclose(dados))

	{

		printf("erro");

	}

}





void writeT (dataFromSOR data)

{



int g,h;

		

				FILE *solucao;

				solucao = fopen ("solucao8x8.txt", "a");

				

				if (solucao == NULL)

				{

				printf ("Houve um erro ao criar o arquivo dados.\n");

				}

			

			fprintf(solucao, "\n\n w = %2.2f \n", data.w);

			

				for (g=2; g<NY+1; g++)

				{

					for (h=2; h<NX+1; h++)

					{

					fprintf(solucao, "%2.7f \t", data.T[g][h]);

					}

				}

					

				if(fclose(solucao)) 

				{

					printf("erro");

				}

}







int main ()

{	



	interface();

	

	

	

	int i, j;

	float w;

	dataFromSOR data;

	FMatrix FM;

	



	double a, b, c, d, e, f; 

							

	a = -K*HY*HY;	

	b = -K*HY*HY;

	c = -K*HX*HX;

	d = -K*HX*HX;

	e = 2*K*(HX*HX + HY*HY);

	f = HX*HX*HY*HY;

	

	FM = findF ();

	

 //	for (w = 0.01 ; w <= 1.99 ; w = w + 0.01)

//	{

	w=1.0;

	data = SOR (w, FM, a, b, c, d, e, f);  

	writeW(data);

	writeT(data);

//	}

		

}	 

removido
(usa Nenhuma)

Enviado em 28/10/2016 - 18:31h

Apenas para não haver dúvidas, o que seria SOR?
É que tem o método das sobre/sub-relaxações sucessivas (SOR/SUR).

----------------------------------------------------------------------------------------------------------------
Nem direita, nem esquerda. Quando se trata de corrupção o Brasil é ambidestro.
(anônimo)