[Lapackpp-cvs] CVS: lapackpp/testing tCmplxSolve.cc,NONE,1.1 .cvsignore,1.1,1.2 Makefile.am,1.4,1.5 lapack++_test,1.3,1.4

[Christian S. <cs...@us...>]

Update of /cvsroot/lapackpp/lapackpp/testing
In directory sc8-pr-cvs1.sourceforge.net:/tmp/cvs-serv23186/testing

Modified Files:
	.cvsignore Makefile.am lapack++_test 
Added Files:
	tCmplxSolve.cc 
Log Message:
2004-03-01  Christian Stimming <sti...@tu...>

        * testing/tCmplxSolve.cc: Add test program for solving
        complex-valued systems.



--- NEW FILE: tCmplxSolve.cc ---
//
//              LAPACK++ 1.1 Linear Algebra Package 1.1
//               University of Tennessee, Knoxvilee, TN.
//            Oak Ridge National Laboratory, Oak Ridge, TN.
//        Authors: J. J. Dongarra, E. Greaser, R. Pozo, D. Walker
//                 (C) 1992-1996 All Rights Reserved
//
//                             NOTICE
//
// Permission to use, copy, modify, and distribute this software and
// its documentation for any purpose and without fee is hereby granted
// provided that the above copyright notice appear in all copies and
// that both the copyright notice and this permission notice appear in
// supporting documentation.
//
// Neither the Institutions (University of Tennessee, and Oak Ridge National
// Laboratory) nor the Authors make any representations about the suitability 
// of this software for any purpose.  This software is provided ``as is'' 
// without express or implied warranty.
//
// LAPACK++ was funded in part by the U.S. Department of Energy, the
// National Science Foundation and the State of Tennessee.

#if 0
#include "lapack++.h"
#endif

#if 1
#include "lafnames.h"       /* macros for LAPACK++ filenames */
#include LA_GEN_MAT_COMPLEX_H
#include LA_VECTOR_DOUBLE_H 
#include "blas++.h"
#include LA_SOLVE_DOUBLE_H
#include LA_GENERATE_MAT_DOUBLE_H
#include LA_EXCEPTION_H
#include LA_UTIL_H
#endif

bool output = true;

double residual(LaGenMatComplex &A, LaVectorComplex &x, 
    const LaVectorComplex& b)
{
    int M = A.size(0);
    int N = A.size(1);

    LaVectorComplex Axb(b);
    Blas_Mat_Vec_Mult(A, x, Axb, 1.0, -1.0); // = A*x-b

    if (output) {
	std::cout << "\tNorm_Inf(A*x-b)" << Norm_Inf(Axb) << std::endl;
	std::cout << "\tNorm_Inf(A) " << Norm_Inf(A) << std::endl;
	std::cout << "\tNorm_Inf(x) " << Norm_Inf(x) << std::endl;
	std::cout << "\tMacheps :" << Mach_eps_double() << std::endl;
    }

    if (M>N)
    {
        LaVectorComplex R(M);

        Blas_Mat_Trans_Vec_Mult(A, Axb, R);
        return Norm_Inf(R) / 
            (Norm_Inf(A)* Norm_Inf(x) * N * Mach_eps_double());

    }
    else
    {
        return Norm_Inf(Axb ) /
                ( Norm_Inf(A)* Norm_Inf(x) * N * Mach_eps_double());
    }
}


int TestGenLinearSolve(int M,int N)
{
    LaGenMatComplex A(M,N);
    LaVectorComplex x(N), b(M);

#ifndef HPPA
    const char fname[] = "TestGenLinearSolve(LaGenMat, x, b) ";
#else
    char *fname = NULL;
#endif

    //char e = 'e';
    double norm;
    double res;

    LaRandUniform(A, -1, 1);

    // save a snapshot of what A looked like before the solution
    LaGenMatComplex old_A = A;



    b = 1.1;

    if (output)
    std::cout << fname << ": testing LaLinearSolve(Gen,...) M= "<< M
        << " N = " << N << std::endl;

    LaLinearSolve(A, x, b);
    LaGenMatComplex diff_A(old_A);
    

    if ( (norm = Norm_Inf( old_A - A)) >  0.0)  // this is an exact test, not
                                         // necessary to worry about
                                         // round-off issues.  We
                                         // are testing to see A was
                                         // overwritten.
    {
        std::cerr << fname << ": overwrote 1st arg.\n";
        std::cerr << "       error norm: " << norm << std::endl;
        exit(1);
    }

    res = residual(A,x,b);
    if (res > 1)
    {
        std::cerr << fname << "resdiual " << res << " is to too high.\n";
        exit(1);
    }

    if (output)
    std::cout << fname << ": LaLinearSolve() success.\n\n";


    // now try the in-place solver


    if (output)
    std::cout << fname << ": testing LaLinearSolveIP(Gen,...) \n";
    LaLinearSolveIP(A, x, b);


    res = residual(old_A, x, b);

    if (res > 1)
    {
        std::cerr << fname << "resdiual " << res << " is to too high.\n";
        exit(1);
    }

    if (output)
    std::cout << fname << ": LaLinearSolveIP() success.\n\n";
    return 0;
}


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

    std::cout.precision(4);
    std::cout.setf(std::ios::scientific, std::ios::floatfield);

    if (argc < 2)
    {
        std::cerr << "Usage " << argv[0] << " N " << std::endl;
        exit(1);
    }
    int N = atoi(argv[1]);
    int M = N;
    if (argc > 2)
	if (std::string(argv[2])=="q")
	    output = false;

    if (output)
	std::cout << "Testing " << N << " x " << N << " system." << std::endl;

    TestGenLinearSolve(M,N);

    return 0;
}


Index: .cvsignore
===================================================================
RCS file: /cvsroot/lapackpp/lapackpp/testing/.cvsignore,v
retrieving revision 1.1
retrieving revision 1.2
diff -C2 -d --brief -r1.1 -r1.2
Files /tmp/cvswz7lnG and /tmp/cvsWDFSbC differ

Index: Makefile.am
===================================================================
RCS file: /cvsroot/lapackpp/lapackpp/testing/Makefile.am,v
retrieving revision 1.4
retrieving revision 1.5
diff -C2 -d --brief -r1.4 -r1.5
Files /tmp/cvsPK23hJ and /tmp/cvs9zCu8E differ

Index: lapack++_test
===================================================================
RCS file: /cvsroot/lapackpp/lapackpp/testing/lapack++_test,v
retrieving revision 1.3
retrieving revision 1.4
diff -C2 -d --brief -r1.3 -r1.4
Files /tmp/cvsQGQJdO and /tmp/cvssx3n7J differ