成都网站开发哪个好,seo网络推广公司报价,网站成立时间,哈尔滨通用建设工程有限公司文章目录前言Matrix MarketMatlab IORead dataWrite data测试C IORead and write dataDownload MatrixIO 代码下载参考网址前言
一般情况集成了一个线性求解器#xff08;Axb#xff09;#xff0c;我们需要验证其性能和精度#xff0c;这时需要大量数据来做验证#xff…
文章目录前言Matrix MarketMatlab IORead dataWrite data测试C IORead and write dataDownload MatrixIO 代码下载参考网址前言
一般情况集成了一个线性求解器Axb我们需要验证其性能和精度这时需要大量数据来做验证 尤其是有不同性质的矩阵 AAA 的数据例如稀疏性 对称性 正定性 对角占优等。
Matrix Market
Matrix Market Matrix Market 网站提供了友好的接口和数据 方便我们验证求解器的精度和性能尤其是对各种性质的矩阵 AAA 都有很多数据 包含不同维数。
Matlab IO
Read data
function [A,rows,cols,entries,rep,field,symm] mmread(filename)
%
% function [A] mmread(filename)
%
% function [A,rows,cols,entries,rep,field,symm] mmread(filename)
%
% Reads the contents of the Matrix Market file filename
% into the matrix A. A will be either sparse or full,
% depending on the Matrix Market format indicated by
% coordinate (coordinate sparse storage), or
% array (dense array storage). The data will be duplicated
% as appropriate if symmetry is indicated in the header.
%
% Optionally, size information about the matrix can be
% obtained by using the return values rows, cols, and
% entries, where entries is the number of nonzero entries
% in the final matrix. Type information can also be retrieved
% using the optional return values rep (representation), field,
% and symm (symmetry).
%mmfile fopen(filename,r);
if ( mmfile -1 )disp(filename);error(File not found);
end;header fgets(mmfile);
if (header -1 )error(Empty file.)
end% NOTE: If using a version of Matlab for which strtok is not
% defined, substitute gettok for strtok in the
% following lines, and download gettok.m from the
% Matrix Market site.
[head0,header] strtok(header); % see note above
[head1,header] strtok(header);
[rep,header] strtok(header);
[field,header] strtok(header);
[symm,header] strtok(header);
head1 lower(head1);
rep lower(rep);
field lower(field);
symm lower(symm);
if ( length(symm) 0 )disp([Not enough words in header line of file ,filename]) disp(Recognized format: )disp(%%MatrixMarket matrix representation field symmetry)error(Check header line.)
end
if ( ~ strcmp(head0,%%MatrixMarket) )error(Not a valid MatrixMarket header.)
end
if ( ~ strcmp(head1,matrix) )disp([This seems to be a MatrixMarket ,head1, file.]);disp(This function only knows how to read MatrixMarket matrix files.);disp( );error( );
end% Read through comments, ignoring themcommentline fgets(mmfile);
while length(commentline) 0 commentline(1) %,commentline fgets(mmfile);
end% Read size information, then branch according to
% sparse or dense formatif ( strcmp(rep,coordinate)) % read matrix given in sparse % coordinate matrix format[sizeinfo,count] sscanf(commentline,%d%d%d);while ( count 0 )commentline fgets(mmfile);if (commentline -1 )error(End-of-file reached before size information was found.)end[sizeinfo,count] sscanf(commentline,%d%d%d);if ( count 0 count ~ 3 )error(Invalid size specification line.)endendrows sizeinfo(1);cols sizeinfo(2);entries sizeinfo(3);if ( strcmp(field,real) ) % real valued entries:[T,count] fscanf(mmfile,%f,3);T [T; fscanf(mmfile,%f)];if ( size(T) ~ 3*entries )message ...str2mat(Data file does not contain expected amount of data.,...Check that number of data lines matches nonzero count.);disp(message);error(Invalid data.);endT reshape(T,3,entries);A sparse(T(:,1), T(:,2), T(:,3), rows , cols);elseif ( strcmp(field,complex)) % complex valued entries:T fscanf(mmfile,%f,4);T [T; fscanf(mmfile,%f)];if ( size(T) ~ 4*entries )message ...str2mat(Data file does not contain expected amount of data.,...Check that number of data lines matches nonzero count.);disp(message);error(Invalid data.);endT reshape(T,4,entries);A sparse(T(:,1), T(:,2), T(:,3) T(:,4)*sqrt(-1), rows , cols);elseif ( strcmp(field,pattern)) % pattern matrix (no values given):T fscanf(mmfile,%f,2);T [T; fscanf(mmfile,%f)];if ( size(T) ~ 2*entries )message ...str2mat(Data file does not contain expected amount of data.,...Check that number of data lines matches nonzero count.);disp(message);error(Invalid data.);endT reshape(T,2,entries);A sparse(T(:,1), T(:,2), ones(entries,1) , rows , cols);endelseif ( strcmp(rep,array) ) % read matrix given in dense % array (column major) format[sizeinfo,count] sscanf(commentline,%d%d);while ( count 0 )commentline fgets(mmfile);if (commentline -1 )error(End-of-file reached before size information was found.)end[sizeinfo,count] sscanf(commentline,%d%d);if ( count 0 count ~ 2 )error(Invalid size specification line.)endendrows sizeinfo(1);cols sizeinfo(2);entries rows*cols;if ( strcmp(field,real) ) % real valued entries:A fscanf(mmfile,%f,1);A [A; fscanf(mmfile,%f)];if ( strcmp(symm,symmetric) | strcmp(symm,hermitian) | strcmp(symm,skew-symmetric) ) for j1:cols-1,currenti j*rows;A [A(1:currenti); zeros(j,1);A(currenti1:length(A))];endelseif ( ~ strcmp(symm,general) )disp(Unrecognized symmetry)disp(symm)disp(Recognized choices:)disp( symmetric)disp( hermitian)disp( skew-symmetric)disp( general)error(Check symmetry specification in header.);endA reshape(A,rows,cols);elseif ( strcmp(field,complex)) % complx valued entries:tmpr fscanf(mmfile,%f,1);tmpi fscanf(mmfile,%f,1);A tmprtmpi*i;for j1:entries-1tmpr fscanf(mmfile,%f,1);tmpi fscanf(mmfile,%f,1);A [A; tmpr tmpi*i];endif ( strcmp(symm,symmetric) | strcmp(symm,hermitian) | strcmp(symm,skew-symmetric) ) for j1:cols-1,currenti j*rows;A [A(1:currenti); zeros(j,1);A(currenti1:length(A))];endelseif ( ~ strcmp(symm,general) )disp(Unrecognized symmetry)disp(symm)disp(Recognized choices:)disp( symmetric)disp( hermitian)disp( skew-symmetric)disp( general)error(Check symmetry specification in header.);endA reshape(A,rows,cols);elseif ( strcmp(field,pattern)) % pattern (makes no sense for dense)disp(Matrix type:,field)error(Pattern matrix type invalid for array storage format.);else % Unknown matrix typedisp(Matrix type:,field)error(Invalid matrix type specification. Check header against MM documentation.);end
end%
% If symmetric, skew-symmetric or Hermitian, duplicate lower
% triangular part and modify entries as appropriate:
%if ( strcmp(symm,symmetric) )A A A. - diag(diag(A));entries nnz(A);
elseif ( strcmp(symm,hermitian) )A A A - diag(diag(A));entries nnz(A);
elseif ( strcmp(symm,skew-symmetric) )A A - A;entries nnz(A);
endfclose(mmfile);
% Done.
Write data
function [ err ] mmwrite(filename,A,comment,field,precision)
%
% Function: mmwrite(filename,A,comment,field,precision)
%
% Writes the sparse or dense matrix A to a Matrix Market (MM)
% formatted file.
%
% Required arguments:
%
% filename - destination file
%
% A - sparse or full matrix
%
% Optional arguments:
%
% comment - matrix of comments to prepend to
% the MM file. To build a comment matrix,
% use str2mat. For example:
%
% comment str2mat( Comment 1 ,...
% Comment 2,...
% and so on.,...
% to attach a date:,...
% [ ,date]);
% If ommitted, a single line date stamp comment
% will be included.
%
% field - real
% complex
% integer
% pattern
% If ommitted, data will determine type.
%
% precision - number of digits to display for real
% or complex values
% If ommitted, full working precision is used.
%if ( nargin 5) precision 16;
elseif ( nargin 4) precision 16;
elseif ( nargin 3) mattype real; % placeholder, will check after FIND-ing Aprecision 16;
elseif ( nargin 2) comment ;% Check whether there is an imaginary part:mattype real; % placeholder, will check after FIND-ing Aprecision 16;
endmmfile fopen([filename],w);
if ( mmfile -1 )error(Cannot open file for output);
end;[M,N] size(A);%%%%%%%%%%%%% This part for sparse matrices %%%%%%%%%%%%%%%%
if ( issparse(A) )[I,J,V] find(A);if ( sum(abs(imag(nonzeros(V)))) 0 )Vreal 0; else Vreal 1; endif ( ~ strcmp(mattype,pattern) Vreal )mattype real; elseif ( ~ strcmp(mattype,pattern) )mattype complex;end
%
% Determine symmetry:
%if ( M ~ N )symm general;issymm 0;NZ length(V);elseissymm 1;NZ length(V);for i1:NZif ( A(J(i),I(i)) ~ V(i) )issymm 0;break;endendif ( issymm )symm symmetric;ATEMP tril(A);[I,J,V] find(ATEMP);NZ nnz(ATEMP);elseisskew 1;for i1:NZif ( A(J(i),I(i)) ~ - V(i) )isskew 0;break;endendif ( isskew )symm skew-symmetric;ATEMP tril(A);[I,J,V] find(ATEMP);NZ nnz(ATEMP);elseif ( strcmp(mattype,complex) )isherm 1;for i1:NZif ( A(J(i),I(i)) ~ conj(V(i)) )isherm 0;break;endendif ( isherm )symm hermitian;ATEMP tril(A);[I,J,V] find(ATEMP);NZ nnz(ATEMP);else symm general;NZ nnz(A);endelsesymm general;NZ nnz(A);endendend% Sparse coordinate format:rep coordinate;fprintf(mmfile,%%%%MatrixMarket matrix %s %s %s\n,rep,mattype,symm);[MC,NC] size(comment);if ( MC 0 )fprintf(mmfile,%% Generated %s\n,[date]);elsefor i1:MC,fprintf(mmfile,%%%s\n,comment(i,:));endendfprintf(mmfile,%d %d %d\n,M,N,NZ);cplxformat sprintf(%%d %%d %% .%dg %% .%dg\n,precision,precision);realformat sprintf(%%d %%d %% .%dg\n,precision);if ( strcmp(mattype,real) )for i1:NZfprintf(mmfile,realformat,I(i),J(i),V(i));end;elseif ( strcmp(mattype,complex) )for i1:NZfprintf(mmfile,cplxformat,I(i),J(i),real(V(i)),imag(V(i)));end;elseif ( strcmp(mattype,pattern) )for i1:NZfprintf(mmfile,%d %d\n,I(i),J(i));end;else err -1;disp(Unsupported mattype:)mattypeend;%%%%%%%%%%%%% This part for dense matrices %%%%%%%%%%%%%%%%
elseif ( sum(abs(imag(nonzeros(A)))) 0 )Areal 0; else Areal 1; endif ( ~strcmp(mattype,pattern) Areal )mattype real;elseif ( ~strcmp(mattype,pattern) )mattype complex;end
%
% Determine symmetry:
%if ( M ~ N )issymm 0;symm general;elseissymm 1;for j1:N for ij1:Nif (A(i,j) ~ A(j,i) )issymm 0; break; endendif ( ~ issymm ) break; endendif ( issymm )symm symmetric;elseisskew 1;for j1:N for ij1:Nif (A(i,j) ~ - A(j,i) )isskew 0; break; endendif ( ~ isskew ) break; endendif ( isskew )symm skew-symmetric;elseif ( strcmp(mattype,complex) )isherm 1;for j1:N for ij1:Nif (A(i,j) ~ conj(A(j,i)) )isherm 0; break; endendif ( ~ isherm ) break; endendif ( isherm )symm hermitian;else symm general;endelsesymm general;endendend% Dense array format:rep array;[MC,NC] size(comment);fprintf(mmfile,%%%%MatrixMarket matrix %s %s %s\n,rep,mattype,symm);for i1:MC,fprintf(mmfile,%%%s\n,comment(i,:));end;fprintf(mmfile,%d %d\n,M,N);cplxformat sprintf(%% .%dg %% .%dg\n, precision,precision);realformat sprintf(%% .%dg\n, precision);if ( ~ strcmp(symm,general) )rowloop j;else rowloop 1;endif ( strcmp(mattype,real) )for j1:Nfor ieval(rowloop):Mfprintf(mmfile,realformat,A(i,j));endendelseif ( strcmp(mattype,complex) )for j1:Nfor ieval(rowloop):Mfprintf(mmfile,cplxformat,real(A(i,j)),imag(A(i,j)));endendelseif ( strcmp(mattype,pattern) )err -2disp(Pattern type inconsistant with dense matrix)elseerr -2disp(Unknown matrix type:)mattypeend
endfclose(mmfile);
测试
% 读取 bcsstk14.mtx 的矩阵信息存储于矩阵 A 中 rows cols 分别是行列
[A,rows,cols,entries,rep,field,symm] mmread(bcsstk14.mtx);
% 将矩阵 A 写到 data.txt 文件中
mmwrite(data.txt, A, field)
[A2,rows,cols,entries,rep,field,symm] mmread(data.txt);
error A - A2;
error_norm norm(error, 2)
C IO
在 网址 C IO 下载 mmio.h 和 mmio.c 文件
Read and write data
#include iostream
#include fstream
#include cstdio
#include cstdlib
#include string.h
#include stdlib.h
#include string
#include map
#include mmio.h
#include stdio.husing namespace std;int main(int argc, char *argv[])
{
#if 1// read datastring fileName ../data/bcsstk14.mtx; int ret_code;MM_typecode matcode;FILE *f;int M, N, nz; int i, *I, *J;double *val;mapunsigned int, mapunsigned int, double A; // Storage in triplet mode row col val ; As id from zeroif ((f fopen(fileName.c_str(), r)) NULL) {cerrread file error, please check.endl;exit(1);}if (mm_read_banner(f, matcode) ! 0){printf(Could not process Matrix Market banner.\n);exit(1);}/* This is how one can screen matrix types if their application *//* only supports a subset of the Matrix Market data types. */if (mm_is_complex(matcode) mm_is_matrix(matcode) mm_is_sparse(matcode) ) {printf(Sorry, this application does not support );printf(Market Market type: [%s]\n, mm_typecode_to_str(matcode));exit(1);}/* find out size of sparse matrix .... */if ((ret_code mm_read_mtx_crd_size(f, M, N, nz)) !0) {exit(1);}/* reseve memory for matrices */I (int *) malloc(nz * sizeof(int));J (int *) malloc(nz * sizeof(int));val (double *) malloc(nz * sizeof(double));/* NOTE: when reading in doubles, ANSI C requires the use of the l *//* specifier as in %lg, %lf, %le, otherwise errors will occur *//* (ANSI C X3.159-1989, Sec. 4.9.6.2, p. 136 lines 13-15) */for (i0; inz; i) {fscanf(f, %d %d %lg\n, I[i], J[i], val[i]);I[i]--; /* adjust from 1-based to 0-based */J[i]--;}if (f ! stdin) {fclose(f);}/************************//* now write out matrix *//************************/mm_write_banner(stdout, matcode);mm_write_mtx_crd_size(stdout, M, N, nz);for (i0; inz; i) {//fprintf(stdout, %d %d %20.19g\n, I[i]1, J[i]1, val[i]);A[I[i]][J[i]] val[i];}#endif#if 0// write dataconst int nz 4;const int M 10;const int N 10;MM_typecode matcode; int I[nz] { 0, 4, 2, 8 };int J[nz] { 3, 8, 7, 5 };double val[nz] {1.1, 2.2, 3.2, 4.4};int i;mm_initialize_typecode(matcode);mm_set_matrix(matcode);mm_set_coordinate(matcode);mm_set_real(matcode);mm_write_banner(stdout, matcode); mm_write_mtx_crd_size(stdout, M, N, nz);/* NOTE: matrix market files use 1-based indices, i.e. first elementof a vector has index 1, not 0. */for (i0; inz; i) {fprintf(stdout, %d %d %10.3g\n, I[i]1, J[i]1, val[i]);}
#endifreturn 0;
}
Download Matrix
在网址的主页面可以看到 Search by matrix properties 在之后的红色框中选择所需要的矩阵性质 回车即可。
最后利用前面介绍的 IO 的代码 即可把矩阵导入到自己的线性求解器中 作为验证。
IO 代码下载
以上的所有代码包括一个测试的例子均可以到我的 GitHub 网页下载其中包含 MakeFile文件。
参考网址
Matrix Market https://math.nist.gov/MatrixMarket/
