Gauss-seidel and ASCII

[Brandybuck]

Hey I'm in mechanical engineering and hate coding lol

if someone could help me with this assignment I'd really appreciate it.

The Program

(1) Create an ASCII data file of the coefficient and right-hand side constants of the form:

15 -3 -1 3800
-3 18 -6 1200
-4 -3 12 2350

(2) Write a function that will implement the Gauss-Seidel algorithm.

(3) Write a main program that will read the data from the ASCII file you created, and then call your Gauss-Seidel Function. After the function is called it will write the data to an ASCII file name results.dat in the following format:

C[1] = 326.2458
C[2] = 242.8018
C[3] = 365.2824

These are the answers you should get if your program is working properly.

(4) Put the program and function in one file.

(5) The GA will use g++ to compile your program, so be sure it compiles with this command.

[lastchance]

I'm guessing that, to within the letters you choose to use for them, your matrix M is

15  -3  -1 
-3  18  -6 
-4  -1  12 

and that your RHS vector b is

3800
1200
2350

so that your matrix equation is

MC = b

So your equations FOR THESE SPECIFIC NUMBERS (and please read the comment about the start of array indices further down this post) are

15 C[1]  - 3 C[2] -    C[3] = 3800
-3 C[1] + 18 C[2] -  6 C[3] = 1200
-4 C[1] -    C[2] + 12 C[3] = 2350

Rearrange to make each diagonal element the subject:

C[1] = ( 3800 + 3 C[2] +   C[3] ) / 15
C[2] = ( 1200 + 3 C[1] + 6 C[3] ) / 18
C[3] = ( 2350 + 4 C[1] +   C[2] ) / 12

Gauss-Seidel just keeps looping through these equations until the maximum change in a loop is less than some small tolerance. Looking at how diagonally-dominant your matrix is, it should converge in no time.


These numbers will help you see what is going on, but you should code it in such a way that:
(1) You can use any numbers in the matrix or RHS (at least those for which G-S converges)
(2) Ideally you should have any number of elements (not just 3)

As an additional warning, you should also be aware that array indices in C++ start from 0, not 1.

If you show some code there are lots of people on this forum who will help you to improve / correct it, but the onus is on you to write some code first.

[Brandybuck]

Hey I wrote some code. currently I've hard coded the array into the program, and I'm testing it by just printing C[0].

C[0] is returning as NaN

#include "stdafx.h"
#include<iostream>
#include<cmath>

using namespace std;
int main()
{
double C[3], er, Cold;
double a[3][4] = { 15, -3, -1, 3800
- 3, 18, -6, 1200
- 4, -3, 12, 2350 };

er = 0.01;
C[0] = a[0][3] / a[0][0];
C[1] = a[1][3] / a[1][1];
C[2] = a[2][3] / a[2][2];
Cold = 326.2458;
while (true)
{
if ((Cold - C[0]) / Cold >= er)
{
C[0] = (a[0][3] - C[1] * a[0][1] - C[2] * a[0][2]) / a[0][0];
C[1] = (a[1][3] - C[0] * a[1][0] - C[2] * a[1][2]) / a[1][1];
C[2] = (a[2][3] - C[2] * a[2][0] - C[1] * a[2][1]) / a[2][2];
}
else {

cout << C[0];
break;

}
}

}

[Brandybuck]

nvm found the error
it was greater than I expected ;)

[Brandybuck]

nvm I had it right the first time still coming out as NaN

[jonnin]

I don't see any way to drive it to NAN. Division by zero and uninitialized values or out of bounds on an aray are the suspects. The code does not change the divisors so divide by 0 seems unlikely. I don't see anything.

you probably want fabs(...) >= er, not 100% sure? Could you be going to underflow via this?

just do while off the >= er. don't need to wrap that in an infinite loop and break.

try looping just a few times, maybe replace while true with for... 100 or something to debug it? How many iterations drive it to NAN?

[Brandybuck]

I rewrote the code to do just 1 iteration without a loop and I'm still getting NaN

// ConsoleApplication1.cpp : Defines the entry point for the console application.
//

#include "stdafx.h"
#include<iostream>
#include<cmath>
#include<thread>

using namespace std;
int main()
{
double C[3];
double a[3][4] = { 15, -3, -1, 3800
- 3, 18, -6, 1200
- 4, -3, 12, 2350 };
C[0] = a[0][3] / a[0][0];
C[1] = a[1][3] / a[1][1];
C[2] = a[2][3] / a[2][2];
C[0] = (((a[0][3]) - (C[1] * a[0][1] - C[2] * a[0][2])) / a[0][0]);
C[1] = (((a[1][3]) - (C[0] * a[1][0] - C[2] * a[1][2])) / a[1][1]);
C[2] = (((a[2][3]) - (C[0] * a[2][0] - C[1] * a[2][1])) / a[2][2]);
cout << C[0];
cout << C[1];
cout << C[2];
}

[jonnin]

your a initialize is all wrong
i think its

double a[3][4] = { {15, -3, -1, 3800},
{- 3, 18, -6, 1200},
{- 4, -3, 12, 2350} };

[Brandybuck]

no I already tested that part by itself. my original C[0] value is correct.

[Brandybuck]

I just tried your suggestion it made an infinite loop.

// ConsoleApplication1.cpp : Defines the entry point for the console application.
//

#include "stdafx.h"
#include<iostream>
#include<cmath>
#include<thread>

using namespace std;
int main()
{
double C[3], er, Cold, i;
double a[3][4] = { {15, -3, -1, 3800},
{-3, 18, -6, 1200},
{-4, -3, 12, 2350} };

er = .1;
C[0] = a[0][3] / a[0][0];
C[1] = a[1][3] / a[1][1];
C[2] = a[2][3] / a[2][2];
Cold = 326.2458;
i = 1;
while (true)
{
if (((Cold - C[0]) / Cold) >= er)
{
C[0] = ((a[0][3] - (C[1] * a[0][1] - C[2] * a[0][2])) / a[0][0]);
C[1] = ((a[1][3] - (C[0] * a[1][0] - C[2] * a[1][2])) / a[1][1]);
C[2] = ((a[2][3] - (C[0] * a[2][0] - C[1] * a[2][1])) / a[2][2]);
cout << i;
i = i + 1;
std::this_thread::sleep_for(1s);
}
else {

cout << C[0];
break;

}
}

}

[lastchance]

Please put your code in code tags so that we can read it properly and test it. I agree with Jonnin: your original initialisation was wrong. You would also be less likely to make mistakes if you did the calculations with proper indexed loops and the code would be easier to read if it wasn't so heavily bracketed.

I don't think it's a good idea to put your RHS vector into a[][]. It makes little mathematical sense.

[lastchance]

Your convergence criteria isn't sensible as it is premised on knowing the solution you are trying to find, and also ignores whether the relative error is positive or negative. Debug your code by doing a fixed number of loops first. There is absolutely no need to use threads or sleep functions. You could simply output all components of C[] at the end of each loop and check by hand what is happening.

[Brandybuck]

#include "stdafx.h"
#include<iostream>
#include<cmath>
#include<thread>

using namespace std;
int main()
{
	double C[3], er, Cold, i;
	double a[3][4] = { {15, -3, -1, 3800},
	{-3, 18, -6, 1200},
	{-4, -3, 12, 2350} };

	er = .1;
	C[0] = a[0][3] / a[0][0];
	C[1] = a[1][3] / a[1][1];
	C[2] = a[2][3] / a[2][2];
	Cold = 326.2458;
	i = 1;
	while (true)
	{
		if (((Cold - C[0]) / Cold) >= er)
		{
			C[0] = ((a[0][3] - (C[1] * a[0][1] - C[2] * a[0][2])) / a[0][0]);
			C[1] = ((a[1][3] - (C[0] * a[1][0] - C[2] * a[1][2])) / a[1][1]);
			C[2] = ((a[2][3] - (C[0] * a[2][0] - C[1] * a[2][1])) / a[2][2]);
			cout << i;
			i = i + 1;
			std::this_thread::sleep_for(1s);
		}
		else {

			cout << C[0];
			break;

		}
	}

}

[Brandybuck]

okay so this code is now working, but I cant get the loop to not be infinite

clear, clc
// ConsoleApplication1.cpp : Defines the entry point for the console application.
//

#include "stdafx.h"
#include<iostream>
#include<cmath>
#include<thread>

using namespace std;
int main()
{
	double C[3];
	double a[3][4] = { { 15, -3, -1, 3800 },
	{ -3, 18, -6, 1200 },
	{ -4, -3, 12, 2350 } };
	C[0] = a[0][3] / a[0][0];
	C[1] = a[1][3] / a[1][1];
	C[2] = a[2][3] / a[2][2];
	
	C[0] = ((a[0][3] - (C[1] * a[0][1] - C[2] * a[0][2])) / a[0][0]);
	C[1] = ((a[1][3] - (C[0] * a[1][0] - C[2] * a[1][2])) / a[1][1]);
	C[2] = ((a[2][3] - (C[0] * a[2][0] - C[1] * a[2][1])) / a[2][2]);
	
	cout << C[0];
	cout << C[1];
	cout << C[2];

[jonnin]

did you try the fabs suggestion?

[Brandybuck]

I'll try it now. didnt see that, thanks

[lastchance]

Well, for now do a fixed number of loops - you can sort out a convergence criterion once you have properly indexed arrays. Here's some suggestions:

- make a[3][3] and have a separate rhs vector, e.g. rhs[3]
- use indexed loops and don't have hard-coded separate calculations for C[0], C[1], C[2]
- your calculations are actually wrong because of overall sign errors caused by your over-bracketing - lines 21, 22, 23 in your last code (which isn't complete BTW)
- a convergence criterion should be based on the change over one iteration - not a presumed solution (which you won't get because the calculation formulae have errors in); there are various ways of getting the 'change' - max change over all components of C[], sum of changes for all components of C[], rms change etc - take your pick

[Brandybuck]

okay so my code seems to be working but its not converging at the expected values. maybe my professor was wrong he makes small mistakes a lot

#include "stdafx.h"
#include<iostream>
#include<cmath>

using namespace std;
int main()
{
	double C[3], cL, i;
	double a[3][4] = { { 15, -3, -1, 3800 },
	{ -3, 18, -6, 1200 },
	{ -4, -3, 12, 2350 } };
	C[0] = a[0][3] / a[0][0];
	C[1] = a[1][3] / a[1][1];
	C[2] = a[2][3] / a[2][2];
	i = 1;
	cL = 2000;
	while (true)
	{
		if (i <= cL)
		{
			C[0] = (a[0][3] - (C[1] * a[0][1] - C[2] * a[0][2])) / a[0][0];
			C[1] = (a[1][3] - (C[0] * a[1][0] - C[2] * a[1][2])) / a[1][1];
			C[2] = (a[2][3] - (C[0] * a[2][0] - C[1] * a[2][1])) / a[2][2];
			i = i + 1;
		}
		else
		{
			cout << C[0];
			cout << C[1];
			cout << C[2];
			break;
		}
	}
	
}

[Brandybuck]

I just checked the answers and he was correct. so I must have a mistake somewhere

[Brandybuck]

hahaha took out some brackets and it runs fine now!! thanks guys I really appreciate your help