//
// 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.
//#include "lapackpp.h"
#define LA_BOUNDS_CHECK
#include <stdlib.h> /* for atoi() and exit() */
#include "lafnames.h"
#include LA_EXCEPTION_H
#include LA_GEN_MAT_DOUBLE_H
#include "laexcp.h"
#include LA_TEMPLATES_H
#include <string>
bool output = true;
int test_scale(int N)
{
LaGenMatDouble bla1, bla2;
la::ones(bla1, N);
la::ones(bla2, N);
bla1 *= 2.0;
bla2.add(1.0);
if (bla1.equal_to(bla2))
if (output)
std::cout << "la::equal is true (correct)" << std::endl;
else {
if (output)
std::cout << "la::equal is false (oops, wrong)" << std::endl;
return 1;
}
return 0;
}
int test_templates()
{
LaGenMatDouble bla, eye;
la::ones(bla, 4, 5);
la::zeros(bla, 3, 3);
la::eye(eye, 3, 3);
LaGenMatDouble foo(la::ones<LaGenMatDouble>(1, 3));
la::from_diag(bla, foo);
if (la::equal(eye, bla))
if (output)
std::cout << "la::equal is true (correct)" << std::endl;
else {
if (output)
std::cout << "la::equal is false (oops, wrong)" << std::endl;
return 1;
}
LaGenMatDouble eye2(la::from_diag(la::ones<LaGenMatDouble>(1, 4)));
if (la::equal(eye2, la::eye<LaGenMatDouble>(4, 4)))
if (output)
std::cout << "la::equal is true (correct)" << std::endl;
else {
if (output)
std::cout << "la::equal is false (oops, wrong)" << std::endl;
return 1;
}
la::trace(bla);
la::rand(bla);
bla = la::linspace<LaGenMatDouble>(1, 49, 10);
LaGenMatDouble blub = la::repmat(bla, 2, 4);
if (output)
std::cout << "repmat: " << blub << "end of blub." << std::endl;
if (output)
std::cout << "Row index 1 of blub: " << blub.row(1) << std::endl;
blub.row(1) = 2.0;
if (output)
std::cout << "Row 1 assigned to 2; blub is: " << blub << std::endl;
if (blub(1,1) != 2.0)
return 1;
return 0;
}
int test_index()
{
LaGenMatDouble A(10,10), B, C;
LaIndex II(1,9,2);
LaIndex JJ(1,1,1);
// Note: The elements of A are not yet initialized; now
// initialize them to zero
A = 0.0;
B.ref(A(II,II)); // B references the following elements of A:
// (1, 1), (1, 3), (1, 5), (1, 7), (1, 9)
// (3, 1), (3, 3), (3, 5), (3, 7), (3, 9)
// (5, 1), (5, 3), (5, 5), (5, 7), (5, 9)
// (7, 1), (7, 3), (7, 5), (7, 7), (7, 9)
// (9, 1), (9, 3), (9, 5), (9, 7), (9, 9)
B(2,3) = 3.1; // Set A(5,7) = 3.1
if (A(5, 7) != 3.1)
return 1;
C.ref(B(LaIndex(2,2,4), JJ));
// C references the following elements of B:
// (2, 1)
// C references the following elements of A:
// (5, 3)
C = 1.1; // Set B(2,1) = A(5,3) = 1.1
if (A(5, 3) != 1.1)
return 1;
if (B(2, 1) != 1.1)
return 1;
if (output) {
std::cout << "Test the index operator; A: the following matrix should be almost zeros" << std::endl;
std::cout << A << std::endl;
std::cout << "Test the index operator; B: This should be a 5x5 matrix" << std::endl;
std::cout << B << std::endl;
std::cout << "Test the index operator; C: This should return '1.1'" << std::endl;
std::cout << C;
C.Info(std::cout) << std::endl;
}
return 0;
}
int test_inject()
{
LaGenMatDouble A(3,3), B(3,3), C(2,2);
A = 0;
// element assignments to submatrix views will work:
A(LaIndex(0,2),LaIndex(0,2))=1;
B(LaIndex(0,2),LaIndex(0,2))=2;
// note: parts of B are still uninitialized
C(LaIndex(0,1),LaIndex(0,1))=3;
A(LaIndex(0,1),LaIndex(0,1)).inject(C);
// B(LaIndex(0,1),LaIndex(0,1))=C;
// ^^ will break due to copy(), but inject() would work.
// use a submatrix view to modify a column of A
LaGenMatDouble D(A(LaIndex(), LaIndex(1)));
D(0,0) = 4;
D(1,0) = 5;
D(2,0) = 6;
// and check whether the resulting matrix is correct
double expected_res[] = { 3, 4, 1,
3, 5, 1,
1, 6, 1};
LaGenMatDouble expected(expected_res, 3, 3, true); // row-ordered
if (output)
std::cout<< " Matrix A: " << A <<std::endl;
//std::cout<<B<<std::endl;
if (!expected.equal_to(A))
return 1;
return 0;
}
int main(int argc, char *argv[])
{
int M, N;
int i,j;
LaException::enablePrint(true);
if (argc <3)
{
std::cout << "Usage " << argv[0] << " M N " << std::endl;
assert(argc >= 3);
exit(1);
}
if (argc > 3)
if (std::string(argv[3])=="q")
output = false;
M = atoi(argv[1]);
N = atoi(argv[2]);
double v[100]; // = {4.0};
for (int k=0;k<100;k++) v[k] = 4.0;
// Test constructors
//
LaGenMatDouble A;
if (output)
std::cout << std::endl << "null consturctor " << std::endl;
if (output)
std::cout << "A(): " << A.info() << std::endl;
if (output)
std::cout << std::endl << "(int, int) constructor " << std::endl;
LaGenMatDouble C(M,N);
if (output)
std::cout << "C(M,N): " << C.info() << std::endl;
// C.debug(1);
if (output)
std::cout << std::endl << "X(const &X) constructor " << std::endl;
LaGenMatDouble D(C); // D is also N,N
if (output)
std::cout << "D(C) :" << D.info() << std::endl;
if (output)
std::cout << std::endl;
LaGenMatDouble K, O(100,100);
LaGenMatDouble L(O(LaIndex(2,4),LaIndex(2,8)));
if (output)
std::cout << "L(O(LaIndex(2,4),LaIndex(2,8)))\n";
L = 0.0;
if (output)
std::cout << "L: " << L.info() << std::endl;
if (output)
std::cout << std::endl <<"K.copy(L) " << std::endl;
K.copy(L);
if (output)
std::cout << "K: " << K.info() << std::endl;
if (output)
std::cout << "K:\n" << K << std::endl;
LaIndex I(std::min(2,M-1),M-1), J(1,N-1); // evens, odd
if (output)
std::cout << std::endl << "create indices I=" << I << ", J=" << J << std::endl;
LaGenMatDouble E(C(I, J));
if (output)
std::cout << std::endl << "X(const &X) constructor with submatrices " << std::endl;
if (output)
std::cout << "E(C(I,J)): " << E.info() << std::endl;
for (j=0;j<N; j++)
for (i=0; i<M; i++)
C(i,j) = i + j/100.0;
if (output) {
std::cout << std::endl;
std::cout << "test operator(int, int)" << std::endl;
std::cout << "Initalize C(i,j) = i + j/100.0 " << std::endl;
std::cout << "C: " << C.info() << std::endl;
std::cout << C << std::endl;
std::cout << std::endl;
std::cout << "operator(LaIndex, LaIndex)" << std::endl;
std::cout << "C(I,J) " << C(I,J).info() << std::endl;
std::cout << C(I,J) << std::endl;
std::cout << std::endl;
std::cout << "test missing indices (default to whole row or column" << std::endl;
std::cout << "C(LaIndex(),J) " << C(LaIndex(),J).info() << std::endl;
std::cout << C(LaIndex(),J) << std::endl;
std::cout << std::endl;
std::cout << "C(I,LaIndex()) " << C(I,LaIndex()).info() << std::endl;
std::cout << C(I,LaIndex()) << std::endl;
std::cout << std::endl;
}
LaGenMatDouble F;
if (output)
std::cout << "F.ref(C(I,J))\n";
F.ref(C(I,J));
if (output)
std::cout << F.info() << std::endl;
F = 4.44;
if (output)
std::cout <<"F:\n" << std::endl;
if (output)
std::cout << F << std::endl;
E.inject(F); // changed due to changed operator= semantics
if (output) {
std::cout << std::endl;
std::cout << "operator=() " << std::endl;
std::cout << "E = F : " << E.info() << C << std::endl;
}
D = C;
if (output) {
std::cout << std::endl;
std::cout << "operator=(const Matrix&) "<< std::endl;
std::cout << "D = C : " << D.info() << std::endl;
std::cout << D << std::endl;
std::cout << std::endl;
std::cout << "test automatic destructuion of temporaries: " << std::endl;
}
LaGenMatDouble B;
for (int c=0; c<10; c++)
{
B.ref(C(I,J));
if (output)
std::cout << "B.ref(C(I,J)): " << B.info() << std::endl;
}
C.ref(C);
if (output) {
std::cout << std::endl;
std::cout <<"test C.ref(C) case works correctly." << std::endl;
std::cout << "C.ref(C) " << C.info() << std::endl;
}
return test_index() +
test_inject() +
test_templates() +
test_scale(N);
}
syntax highlighted by Code2HTML, v. 0.9.1