# Code: One dimensional steady state conduction with heat generation

(Difference between revisions)
 Revision as of 08:05, 9 June 2006 (view source)Praveen (Talk | contribs)← Older edit Latest revision as of 12:32, 19 December 2008 (view source)Peter (Talk | contribs) m (Reverted edits by DommoNrela (Talk) to last version by Praveen) (One intermediate revision not shown)

## Latest revision as of 12:32, 19 December 2008

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

#define MAXROWS 10000

double **GetMatrix (double, double, double, double, double, int);
void squareoutput (double **, int);
void nonsquareoutput (double **, int, int);
double **triangularize (double *[MAXROWS], int, int);
double *returnsolvector (double **c, int nrows);
void writeoutputvector (double *p, int nrows);

int
main (int argc, char **argv)
{

/*
* Variables chosen:
* GAMMA:  Thermal Conductivity
* L    : Length of the domain
* Dx   : Node spacing
* */

double GAMMA, L, Dx, Ta, Tb, Area;
int N, row;
double **temp;
double **aug;
double *solvector;

printf ("\nThis program works for only simple 1-D conduction");
printf ("Is thermal Conductivity constant?");

printf ("\nEnter the Value for Heat conductivity");
scanf ("%lf", &GAMMA);
printf ("\nEnter the Length of the computational domain");
scanf ("%lf", &L);
printf ("\nEnter the cross-sectional area of the computational domain");
scanf ("%lf", &Area);
printf ("\nEnter the temperature at the ends of the domain");
scanf ("%lf %lf", &Ta, &Tb);
printf ("\nEnter the number of Nodes u want to be chosen for analysis");
scanf ("%d", &N);

temp = (double **) malloc (N * sizeof (double *));
for (row = 0; row < N; row++)
temp[row] = (double *) malloc (N * sizeof (double));

aug = (double **) malloc ((N + 1) * sizeof (double));
for (row = 0; row < N; row++)
aug[row] = (double *) malloc ((N + 1) * sizeof (double));

solvector = (double *) malloc (N * sizeof (double));

Dx = L / (N);

printf ("\nThe Node spacing is: %lf", Dx);

aug = GetMatrix (GAMMA, Dx, Ta, Tb, Area, N);
aug = triangularize (aug, N, N + 1);
nonsquareoutput (aug, N, N + 1);
printf ("\n\n");
solvector = returnsolvector (aug, N);
writeoutputvector (solvector, N);
printf ("\n\n\n");
return 0;
}

double **
GetMatrix (double GAMMA, double Dx, double Ta, double Tb, double Area, int N)
{
int row, col;
double DxW, DxE;
double Aw, Ae, Ap;

/* Memory Allocation for the matrix */

double **matrix;
double *solvector;
double **aug;

solvector = (double *) malloc (N * sizeof (double));

matrix = (double **) malloc (N * sizeof (double));
for (row = 0; row < N; row++)
matrix[row] = (double *) malloc (N * sizeof (double));

aug = (double **) malloc ((N + 1) * sizeof (double));
for (row = 0; row < N; row++)
aug[row] = (double *) malloc ((N + 1) * sizeof (double));

/* Uniform Mesh has been chosen */

DxW = DxE = Dx;

Aw = (GAMMA / DxW) * Area;
Ae = (GAMMA / DxE) * Area;

/* Generation of the Matrix */

row = 0;
for (col = 0; col < N; col++) {
if (col == row)
*(*(matrix + row) + col) = (Ae + 2 * Aw);
else if (col == row + 1)
*(*(matrix + row) + col) = -Ae;
else
*(*(matrix + row) + col) = 0;
}

for (row = 1; row < N - 1; row++) {
for (col = 0; col < N; col++) {
if (col == row + 1)
*(*(matrix + row) + col) = -Aw;
else if (col == row - 1)
*(*(matrix + row) + col) = -Ae;
else if (col == row)
*(*(matrix + row) + col) = (Aw + Ae);
else
*(*(matrix + row) + col) = 0;
}
}

row = N - 1;
for (col = 0; col < N; col++) {
if (col == row)
*(*(matrix + row) + col) = (2 * Ae + Aw);
else if (col == row - 1)
*(*(matrix + row) + col) = -Aw;
else
*(*(matrix + row) + col) = 0;
}

/* Getting the solution vector */

*(solvector + 0) = 2 * Aw * Ta;

for (row = 1; row <= N - 2; row++)
*(solvector + row) = 0;

*(solvector + N - 1) = 2 * Ae * Tb;

/* Getting the augmented matrix */

for (row = 0; row < N; row++) {
for (col = 0; col < N; col++) {
*(*(aug + row) + col) = *(*(matrix + row) + col);
}
}

for (row = 0; row < N; row++)
*(*(aug + row) + N) = *(solvector + row);

return (aug);
}

/* Displaying the matrix */

void
squareoutput (double **c, int nrows)
{
int row, col;
for (row = 0; row < nrows; row++) {
printf ("\n");
printf ("\n");
for (col = 0; col < nrows; col++) {
printf ("%8.2lf", *(c[row] + col));
}
}
return;
}

void
nonsquareoutput (double **c, int nrows, int ncols)
{
int row, col;
for (row = 0; row < nrows; row++) {
printf ("\n");
printf ("\n");
for (col = 0; col < ncols; col++) {
printf ("%12.2lf", *(c[row] + col));
}
}
return;
}

double **
triangularize (double *a[MAXROWS], int nrows, int ncols)
{
int row, col, k, temp = 0;
double coef;
double **tri;
tri = (double **) malloc (ncols * sizeof (double *));
for (row = 0; row < nrows; row++)
tri[row] = (double *) malloc (ncols * sizeof (double));

{
int row, col;
for (row = 0; row < nrows; row++)
for (col = 0; col < ncols; col++)
*(tri[row] + col) = *(a[row] + col);
}

for (temp = 0; temp < nrows - 1; temp++) {
if ((tri[temp] + temp) == 0)
continue;
else {
for (row = temp; row < nrows; row++) {
for (k = row + 1; k < nrows; k++) {
if (*(tri[k] + temp) == 0)
continue;
else
coef = (*(tri[k] + temp)) / (*(tri[temp] + temp));
for (col = 0; col < ncols; col++)
*(tri[k] + col) =
*(tri[k] + col) - (coef * *(tri[row] + col));
}
}
}
}

for (temp = nrows - 1; temp > 0; temp--) {
if ((tri[temp] + temp) == 0)
continue;
else {
for (row = temp; row > 0; row--) {
for (k = row - 1; k >= 0; k--) {
if (*(tri[k] + temp) == 0)
continue;
else
coef = (*(tri[k] + temp)) / (*(tri[temp] + temp));
for (col = 0; col < ncols; col++)
*(tri[k] + col) =
*(tri[k] + col) - (coef * *(tri[row] + col));
}
}
}
}

return (tri);
}

double *
returnsolvector (double **c, int nrows)
{
int temp;
double *solvector;
solvector = (double *) malloc (nrows * sizeof (double));

for (temp = 0; temp < nrows; temp++) {
*(solvector + temp) = *(*(c + temp) + nrows) / *(*(c + temp) + temp);
}
return (solvector);
}

void
writeoutputvector (double *p, int nrows)
{
int temp;
for (temp = 0; temp < nrows; temp++)
printf ("\nTemperature at Node %d = %lf", temp + 1, *(p + temp));
}
```