Method of elimination of Gauss with pivoting partial

Solve Linear Equation in format Ax=b

Description

Solve Linear Equation in format Ax=b with method of elimination of Gauss with pivoting partial. Of the 6 file MyGaussSolve2 is the main. Other are auxiliary function. test.m tests the program on some value.

• MATLAB release: MATLAB 6.0 (R12)

[John D'Errico]

For some reason, this does not download. That is for the best anyway. Waste not your time on this.

[Giovanni Di Cecca, Virginia Bellino]

For documentation in Italian Language

http://www.dicecca.net/giovanni/biblo/calcolo/

[C. S.]

So John, you were unable to download the file? How then can you review it? Are you comfortable with giving a bad review when you have not even read the code much less run it? I just lost all respect for your reviews.

[Duane Hanselman]

John: C.S. may have lost all respect for your reviews, but I found it. The zip file is now downloadable and your statement "Waste not your time on this." fits perfectly. Non standard file structure, absolutely no vectorization, not written in English, multiple files when subfunctions are more appropriate for helper files not meant to be called as separate functions, not needed when \ is more appropriate in all cases, not educational even though that could be the only purpose of this submission. I deleted the files.
C.S., download the files, then post another review.

[C. S.]

Duane, I merely objected to the blind review John offered. I made no comments on the code since I too could not download it. TMW had some problem with code and artwork links to recent submissions this morning.

[John D'Errico]

I was able to download it finally after several (4) attempts today. (My initial assumption was the author had failed to upload a proper file.) So here is a better assessment of this sorry mess.

This is simple gaussian elimination with pivoting. Nothing special. Its fully looped, so it will be terribly slow. Lots of crossed function calls, so there is much additional function overhead. I can't compare it for speed, as it fails to run - see below.

The help is useless, as there is none. You need to edit the function to learn how to use it. I hope you know Italian.

When you do try to run this, as I did, it fails to run:

tic,x=myGaussSolve2(A,b);toc
??? Error: File: myGaussSolve2.m Line: 115 Column: 5
A BREAK statement appeared outside of a loop. Use RETURN instead.

There actually are some positive aspects to this code. It checks for errors in the input. It has frequent comments that explain what is done. It even appears to check for zero pivots. But, it is still no more than Gaussian elimination.

This is nothing worth wasting your time to use. Backslash is a far better choice if you seriously want to solve a system of equations. If this actually ran, I might give it a higher rating because of the few positive things I saw. But code that fails to run and has no useable help gets a "poor" rating.

As for Spam's comment, at least I have the courage to put my name on my reviews.

[C. S.]

At least you read the code. Normally that is done before the review is offered but don't let that stop you.

[John D'Errico]

I did read the code in depth. My review is a fair assessment of the work posted. Did you bother to read the code? And do you have the courage to do anything more than post anonymous insults? Apparently not.

[Marco -]

Very excellent.
This project improve a modular approach.

[Tim Davis]

This code is broken, for all matrices. It fails with the error message "A BREAK statement appeared outside a loop ... use RETURN instead". After replacing all relevant break statements with return statements, the test still fails (test 2) with

Error in ==> myGaussSolve2 at 36
iflag = errore (A,B);

??? Output argument "x" (and maybe others) not assigned during call to
"/amd/netapp3/vol/homes11/davis/Desktop/Metodo_di_eliminazioni_di_Gauss_con_pivoting_parziale/myGaussSolve2.m
(myGaussSolve2)".

Error in ==> test at 69
[A,B,pivot,det,iflag,x]=myGaussSolve2(a,b)

the error handling code is broken here. This is for a matrix which is singular to working precision; the code checks to see if the pivot is zero and if so, it returns without setting x to anything. That causes MATLAB to fail. An "error" statement should appear here instead of returning an error flag.

The code is unnecesarily wordy. For example, the maximus function returns the largest entry in A(pivot(k:n),k). It uses a for loop, and a total of 9 statements. It can be done in 2 statements:

[maxval,ip] = max (abs (A (pivot (k:n), k))) ;
ip = ip + k-1 ;

in the description it needs to be pointed out
that this program is for illustration only, and not for production use (backslash or lu is much better).

The comments are in Italian, which itself is not a problem. However, I can tell that the comments do not meet MATLAB standards. "help myGaussSolve2" should explain the inputs and outputs and give examples of its use. Instead it reports:

   Progetto Metodo di Gauss con Pivoting Parziale
 
              Programma elaborato da
 
       Giovanni DI CECCA & Virginia BELLINO
            50 / 887 408 / 466
 
              http://www.dicecca.net

which even if I could read Italian it doesn't tell me how to use the code.

It has a vector that grows in a loop with this absurd loop in the initialization:

% Crea il vettore colonna pivot
for i=1:n
   
   pivot(i,1)=[i];
   
end

... where the pivot vector is not created at this point. Umm, this is just pivot=1:n. Quite an absurd loop here.

This code is otherwise highly commented. However, it looks like something written first in Fortran and then translated directly into MATLAB. No vector statements. Errors returned with "iflag" output arguments just like (say) a Fortran subroutine.

This function also attempts to compute the determinant, but it makes no accomadation for overflow or underflow (to be fair, neither does the MATLAB det, but det is a pretty useless function anyway).

The argument ordering is bizarre:

[A,B,pivot,det,iflag,x] = myGaussSolve2(A,B)

A is returned as the upper triangular factor.
B is the solution B=L\B where L is the (unsaved) lower triangular factor. pivot is the row permutation, iflag is the error code, and finally x is the solution to Ax=b.

Another coding gem:

pivot = scambio (pivot, ip, k)

which does uses 3 lines of code that could be replaced with the one-liner:

pivot ([ip k]) = pivot ([k ip])

that is, it just swaps pivot(ip) and pivot(k). As I said, this looks just like a Fortran code directly translated into MATLAB, line for line. Rather a bit too pedantic even for educational use.

For a much better LU factorization / Gaussian elimination code, see Cleve Moler's NCM (lutx.m) and Nick Higham's toolbox. Unless this is fixed, skip it: it doesn't work and those other codes are far easier to read.

[Tim Davis]

I should add that I tried it in MATLAB 7.5; perhaps it works on MATLAB R12.

[Giovanni Di Cecca]

Dear Professor,
  

   I've re-tested my code, with test.m file, and it work good.

This code is works on MATLAB Student Version 6, Ver. 6.0.0.42a Release 12.

You use an AMD Processor (probably), I've used a Pentium 4 Processor

Just for note, I've included a Patch for JRE on P4 Processor to run MATLAB on my PC

(http://www.dicecca.net/public/patch/Matlab_R12_sv_JRE_Patch_P4.exe)

I haven't seen MATLAB 7.5

For the code I've read, with high interests the comment that made, and in particular a singular line code

pivot ([ip k]) = pivot ([k ip])

that include my three line code.

Fortran is a language that I know named, but I haven't studied and used.

The documentation is in Italian language, for one reason: it is one of three projects that Virginia and I have made

for Examination of numerical calculation (in italian is called Calcolo Numerico).

The "Rather a bit too pedantic even for educational use" is simply just what we wanted my Prof!

I searched on the Internet, and I am happy that a person with his experience, it gave me good advice.

I will treasure his advice for a possible (though unlikely) new version.

As sayd the Romans: <<Ubi maior, minor cessat!>> :-)

Best Regards

Giovanni Di Cecca

[Tim Davis]

Please don't feel intimidated by my status as a Professor. That's not my intent. I'm just trying to point out a few things that can improve your code.

Including comments in Italian is fine - it's just that even in Italian the comments are insufficient. See for example "help chol"; it gives you a short description of what the function does, and tells you what the inputs and outputs are.

Yes, I'm sure the code works in MATLAB 6.0, but the "break" statements are a problem. They should be replaced with "return". Also, when you detect an error you should call "error", not return an error "iflag", and particular you should not return without assigning all of your output arguments. That does cause MATLAB 7.x to generate a true error.

The AMD vs Pentium question is not the issue. I do not have an AMD processor. It's a MATLAB 6.x vs 7.x issue; the former is more forgiving about these kinds of errors, it seems.

[Edony Princivil]