?? goodtutorialonc++templateforbeginners.txt
字號:
l[i]=a[i]/u[i-1];
u[i]=b[i]-(l[i])*(c[i-1]);
}
DOUBLE* y=new DOUBLE[nRow];
y[0]=pb[0];
for(i=1;i<=nRow-1;i++)
y[i]=pb[i]-(l[i])*(y[i-1]);
//DOUBLE* x=new DOUBLE[nRow];
pb[nRow-1]=y[nRow-1]/u[nRow-1];
for(i=nRow-2;i>=0;i--)
pb[i]=(y[i]-(c[i])*(pb[i+1]))/u[i];
delete [] a;
delete [] b;
delete [] c;
delete [] l;
delete [] u;
delete [] y;
}
BOOL CMatrix::ReadMatrix(char c)
{
//從文件以文本形式讀入數據到數組;
//檢查隔離符是否正確
if(c!=' ' && c!=',' && c!='$' && c!=' ')
{
AfxMessageBox("數據間隔符只能用',','$',空格和Tab鍵!");
exit(0);
}
int i,j;
typedef CArray<DOUBLE,DOUBLE> DblArray;
DblArray p;
static char BASED_CODE szFilter[] = "Ascii Files(*.MAsc)|*.MAsc|";
//static char BASED_CODE szFilter[] = "Chart Files (*.xlc)|*.xlc|Worksheet Files (*.xls)|*.xls|Data Files (*.xlc;*.xls)|*.xlc; *.xls|All Files (*.*)|*.*||";
CFileDialog dlg(true,"MAsc",NULL,OFN_HIDEREADONLY | OFN_OVERWRITEPROMPT,szFilter);
if(dlg.DoModal()==IDOK)
{
CFile file(dlg.GetPathName(),CFile::modeRead);
CArchive ar(&file,CArchive::load);
CString str;
CStringArray strA;
ar.ReadString(str);
while(str!="")
{
strA.Add((LPTSTR)(LPCTSTR)str);
ar.ReadString(str);
}
for(i=0;i<strA.GetSize();i++)
{
CString str=strA[i];
CString strL,strR;
int k,l;
str.TrimLeft();
str.TrimRight();
l=str.GetLength();
k=str.Find(c);
while(k!=-1)
{
strL=str.Left(k);
strR=str.Right(l-k-1);
str=strR;
str.TrimLeft();
str.TrimRight();
//strL.TrimLeft();//
strL.TrimRight();
DOUBLE t=atof((LPTSTR)(LPCTSTR)strL);
p.Add(t);
strL=strR;
strL.TrimLeft();
strL.TrimRight();
l=strL.GetLength();
k=strL.Find(c);
}
DOUBLE t;
strL.TrimLeft();
strL.TrimRight();
t=atof((LPTSTR)(LPCTSTR)strL);
p.Add(t);
}
int m,n,L;
m=strA.GetSize();
L=p.GetSize();
if((L%m)!=0)
{
AfxMessageBox("數據讀取錯誤!");
exit(0);
}
else
n=L/m;
allocMemory(m,n);
for(i=0;i<m;i++)
for(j=0;j<n;j++)
{
Val[i][j]=p[n*i+j];
DOUBLE tt=Val[i][j];
}
return true;//讀取到數據;
}
else
return false;//用戶選擇了取消;
}
DOUBLE CMatrix::Bspline3(DOUBLE xi, DOUBLE Vx[], DOUBLE Vy[], int n)
{
//n為給定的插值節點總數目;
//先判斷計算點xi是否在插值節點范圍內,并找出所屬的區間位置;
int i;
if(xi>Vx[n-1] || xi<Vx[0])
{
AfxMessageBox("要計算的點不在已給插值點之間!");
exit(0);
}
if(fabs(xi-Vx[0])<1e-20)
return Vy[0];
if((fabs(xi-Vx[n-1])<1e-20))
return Vy[n-1];
int j;
for(i=0;i<=n-2;i++)
{
if(xi<Vx[i+1])
{
j=i;
i=n-1;
}
}
//n為給定的插值節點的數目;
DOUBLE* h=new DOUBLE[n-1];
DOUBLE* f=new DOUBLE[n-1];
DOUBLE* u=new DOUBLE[n-1]; //u[0] is useless;
DOUBLE* ff=new DOUBLE[n-1];//ff[0] is useless;
DOUBLE* d=new DOUBLE[n-1];// d[0] is useless;
DOUBLE* r=new DOUBLE[n-1];
/* -------------------------用GetParam()代替--------------------------------------
for(i=0;i<=n-2;i++)
{
h[i]=Vx[i+1]-Vx[i];
f[i]=(Vy[i+1]-Vy[i])/h[i];
}
for(i=1;i<=n-2;i++)
{
u[i]=h[i-1]/(h[i-1]+h[i]);
r[i]=1-u[i];
ff[i]=(f[i-1]-f[i])/(Vx[i-1]-Vx[i+1]);
d[i]=6*ff[i];
}
---------------------------------------------------------------------------*/
GetParam(Vx,Vy,n,h,f,u,r,ff,d);
CMatrix A(n-2,n-2);
DOUBLE* b=new DOUBLE[n-2];
DOUBLE* M=new DOUBLE[n];
/*--------------------用GetM()代替;---------------------------------
A(0,0)=2;A(0,1)=r[1];
for(i=1;i<=n-3;i++)
{
A(i,i)=2;
A(i,i+1)=r[i+1];
A(i,i-1)=u[i+1];
}
A(n-2,n-2)=2;A(n-2,n-3)=u[n-1];
-------------------------------------------------------------------*/
GetA(u,r,d,A);
//自然邊界條件;M0=Mn=0;
b[0]=d[1];
for(i=1;i<=n-4;i++)
b[i]=d[i+1];
b[n-3]=d[n-2];
A.chase_EQ(b);//用追趕法求解代數方程組;
M[0]=0;
for(i=1;i<=n-2;i++)
M[i]=b[i-1];
M[n-1]=0;
////////////////////////////////
DOUBLE t1,t2,t3,t33,t44,t4;
t1=(Vx[j+1]-xi)*(Vx[j+1]-xi)*(Vx[j+1]-xi);
t1=t1*(M[j])/(6*h[j]);
t2=(xi-Vx[j])*(xi-Vx[j])*(xi-Vx[j]);
t2=t2*(M[j+1])/(6*h[j]);
t3=(Vx[j+1]-xi)/h[j];
t33=Vy[j]-(M[j])*(h[j])*(h[j])/6;
t3=t3*t33;
t4=(xi-Vx[j])/(h[j]);
t44=Vy[j+1]-(M[j+1])*(h[j])*(h[j])/6;
t4=t4*t44;
delete [] h;
delete [] f;
delete [] u;
delete [] ff;
delete [] d;
delete [] r;
delete [] b;//
delete [] M;
double t5=t1+t2+t3+t4;
return t5;
}
DOUBLE* CMatrix::Bspline3(DOUBLE xi[], int m, DOUBLE Vx[], DOUBLE Vy[], int n)
{
DOUBLE* s=new DOUBLE[m];
double t;
int i;
for(i=0;i<=m-1;i++)
{
t=Bspline3(xi[i],Vx,Vy,n);
s[i]=t;
}
return s;
}
void CMatrix::allocMemory(int m, int n)
{
//先釋放Val已經分配的內存,然后按照制定的大小重新分配內存;
/*---------------------------------------------------------------------
注意:實際上,不能夠根據Val是否為NULL來判斷Val是否已經被分配過內存;
比如:假如CMatrix類包含一個系統提供的默認的拷貝構造函數
(即不進行任何初始化的操作的構造函數),那么,Val的值將是一個隨機值,
如果對象是用此默認構造函數生成,則該對象的Val!=NULL,但是此時,Val
被沒有被分配過內存;
至于在本程序中,在CMatrix類的所有的構造函數中都對Val進行了初始化,
故可以根據Val是否為空來判斷Val是否已經分配過內存;
----------------------------------------------------------------------*/
//先調用freeMemory()函數;
//freeMemory()函數會根據Val是否為空來檢查Val是否已經被分配了內存;
freeMemory();
nRow=m;
nCol=n;
Val=new DOUBLE* [nRow];
for(int i=0;i<nRow;i++)
Val[i]=new DOUBLE [nCol];
}
BOOL CMatrix::ReadMatrix(CString strFilePath, char c)
{
//從文件以文本形式讀入數據到數組;
//檢查隔離符是否正確
if(c!=' ' && c!=',' && c!='$' && c!=' ')
{
AfxMessageBox("數據間隔符只能用',','$',空格和Tab鍵!");
exit(0);
}
int i,j;
typedef CArray<DOUBLE,DOUBLE> DblArray;
DblArray p;
CFile file(strFilePath,CFile::modeRead);
CArchive ar(&file,CArchive::load);
CString str;
CStringArray strA;
if((ar.ReadString(str))==false)
{
AfxMessageBox("文件讀取出錯!檢查文件名和文件路徑!");
exit(0);
}
while(str!="")
{
strA.Add((LPTSTR)(LPCTSTR)str);
ar.ReadString(str);
}
for(i=0;i<strA.GetSize();i++)
{
CString str=strA[i];
CString strL,strR;
int k,l;
str.TrimLeft();
str.TrimRight();
l=str.GetLength();
k=str.Find(c);
while(k!=-1)
{
strL=str.Left(k);
strR=str.Right(l-k-1);
str=strR;
str.TrimLeft();
str.TrimRight();
//strL.TrimLeft();//
strL.TrimRight();
DOUBLE t=atof((LPTSTR)(LPCTSTR)strL);
p.Add(t);
strL=strR;
strL.TrimLeft();
strL.TrimRight();
l=strL.GetLength();
k=strL.Find(c);
}
DOUBLE t;
strL.TrimLeft();
strL.TrimRight();
t=atof((LPTSTR)(LPCTSTR)strL);
p.Add(t);
}
int m,n,L;
m=strA.GetSize();
L=p.GetSize();
if((L%m)!=0)
{
AfxMessageBox("數據讀取錯誤!");
exit(0);
}
else
n=L/m;
allocMemory(m,n);
for(i=0;i<m;i++)
for(j=0;j<n;j++)
{
Val[i][j]=p[n*i+j];
DOUBLE tt=Val[i][j];
}
return true;//讀取到數據;
}
void CMatrix::GSElimin_EQ(DOUBLE *pb,DOUBLE TOL)
{
/*-------------------------------------------------------
高斯消去法解線性代數方程組;
Ax=b; pb為數組b;;
A:方陣;
x為代求的解,是與b同維數的數組;
返回值保存在參數pb中;
--------------------------------------------------------*/
if(nRow!=nCol)
{
AfxMessageBox("高斯消元法只適用于方陣!");
exit(0);
}
/*---------------------------------開始消元-----------------------*/
int i;
int iM;
for(i=0;i<=nRow-2;i++)
{
iM=findMainElem(i);
if(iM==i && fabs(Val[i][i]-0)<TOL)
{
AfxMessageBox("矩陣為奇異矩陣,無法求得唯一解!");
exit(0);
}
else if(iM!=i)
{
exchangeRow(i,iM);
DOUBLE bt;
bt=pb[i];
pb[i]=pb[iM];
pb[iM]=bt;
}
int m,n;
DOUBLE l;
for(m=i+1;m<=nRow-1;m++)
{
l=-Val[m][i]/Val[i][i];
Val[m][i]=0;
for(n=i+1;n<=nCol-1;n++)
{
Val[m][n]=Val[m][n]+l*Val[i][n];
}
pb[m]=pb[m]+l*pb[i];
}
}
if(fabs(Val[nRow-1][nCol-1])<TOL)
{
AfxMessageBox("矩陣奇異,無法求得唯一解 !");
exit(0);
}
/*-------------------------消元完畢----------------------*/
/********************開始求解********************************/
DOUBLE* x=new DOUBLE[nRow];
x[nRow-1]=pb[nRow-1]/Val[nRow-1][nRow-1];
int kk;
for(kk=nRow-2;kk>=0;kk--)
{
DOUBLE sum=0;
for(int jj=kk+1;jj<=nRow-1;jj++)
{
sum=sum+Val[kk][jj]*x[jj];
}
sum=pb[kk]-sum;
x[kk]=sum/Val[kk][kk];
}
/********************求解完畢********************************/
for(i=0;i<=nRow-1;i++)
pb[i]=x[i];
delete [] x;
}
void CMatrix::exchangeRow(int i, int j)
{
//僅僅用于高斯消元;
int k;
DOUBLE t;
for(k=i;k<=nCol-1;k++)
{
t=Val[i][k];
Val[i][k]=Val[j][k];
Val[j][k]=t;
}
}
int CMatrix::findMainElem(int row)
{
//僅僅用于高斯消元;
int i;
DOUBLE max;
int n;
max=Val[row][row];
n=row;
for(i=row+1;i<=nRow-1;i++)
{
if(fabs(Val[i][row])>fabs(max))
{
max=Val[i][row];
n=i;
}
}
return n;
}
BOOL CMatrix::SaveMatrix(BOOL bAscii)
{
static char BASED_CODE szFilter[]="Ascii Files(*.MAsc)|*.MAsc|Binary Files(*.MBin)|*.MBin|";
CFileDialog dlg(false,"MatrixFile.MAsc",NULL,OFN_HIDEREADONLY | OFN_OVERWRITEPROMPT,szFilter);
if(dlg.DoModal()==IDOK)
{
CString str=dlg.GetPathName();
CFile file(str,CFile::modeCreate|CFile::modeReadWrite);
CArchive ar(&file,CArchive::store);
int i,j;
CString s;
CString s1;
if(bAscii)
{
for(i=0;i<=nRow-1;i++)
{
s="";
for(j=0;j<=nCol-1;j++)
{
s1.Format("%.8lf ",Val[i][j]);
s=s+s1;
}
ar.WriteString(s);
ar.WriteString("\r\n");//回車換行符號;
}
return true;
}
if(!bAscii)
{
Serialize(ar);
return true;
}
}
return false;
}
void CMatrix::GetSize(int &m, int &n)
{
m=nRow;
n=nCol;
}
BOOL CMatrix::SaveMatrix(CString strPath, BOOL bAscii)
{
CFile file(strPath,CFile::modeCreate|CFile::modeReadWrite);
CArchive ar(&file,CArchive::store);
if(bAscii)
{
int i,j;
CString s;
CString s1;
for(i=0;i<=nRow-1;i++)
{
s="";
for(j=0;j<=nCol-1;j++)
{
s1.Format("%.8lf ",Val[i][j]);
s=s+s1;
}
ar.WriteString(s);
ar.WriteString("\r\n");//回車換行符號;
}
return true;
}
else
{
Serialize(ar);
return true;
}
return false;
}
int CMatrix::GetRowSize()
{
return nRow;
}
int CMatrix::GetColSize()
{
return nCol;
}
BOOL CMatrix::LoadMatrix(CString strPath)
{
CFile file(strPath,CFile::modeRead);
CArchive ar(&file,CArchive::load);
Serialize(ar);
return true;
}
BOOL CMatrix::LoadMatrix()
{
static char BASED_CODE szFilter[] = "Ascii Files(*.MAsc)|*.MAsc|Binary Files(*.MBin)|*.MBin|";
CFileDialog dlg(true,NULL,"2進制數據文件.MBin",OFN_HIDEREADONLY | OFN_OVERWRITEPROMPT,szFilter);
if(dlg.DoModal()==IDOK)
{
CString str=dlg.GetPathName();
CFile file(str,CFile::modeRead);
CArchive ar(&file,CArchive::load);
Serialize(ar);
return true;
}
return false;
}
還有很多算法,比如矩陣求逆,矩陣求秩,矩陣分解,子矩陣的提取,高斯-賽德爾迭代法等暫時還沒有添加進去,有空慢慢寫。
Powered by BlogDriver 2.1
?? 快捷鍵說明
復制代碼
Ctrl + C
搜索代碼
Ctrl + F
全屏模式
F11
切換主題
Ctrl + Shift + D
顯示快捷鍵
?
增大字號
Ctrl + =
減小字號
Ctrl + -