1 double f(double& d) 2 { 3 d*=2; 4 return d; 5 } 6 7 void f1(double&d) 8 { 9 d*=2; 10 } 11 12 double f2(double& d) 13 { 14 d*=2; 15 return d; 16 }
习题1
1 class f3{ 2 int i; 3 public: 4 f3(int i1):i(i1){} 5 //double operator()(double d){return d*i;}//三者不能共存 6 //void operator()(double&d){d*=i;} 7 double operator()(double&d){d*=i;return d;} 8 };
习题2
1 #include <iostream> 2 #include <vector> 3 #include <string> 4 #include <list> 5 #include<fstream> 6 #include <set> 7 #include<algorithm> 8 #include<stdexcept> 9 #include <map> 10 #include<Matrix.h> 11 #include<MatrixIO.h> 12 using namespace Numeric_lib; 13 using namespace std; 14 15 typedef Matrix<double,2>Matrix1; 16 typedef Matrix<double,1>vector1; 17 18 void elim_with_partial_pivot(Matrix1&A,vector1&b) 19 { 20 const Index n=A.dim1(); 21 22 for(Index j=0;j<n;++j) 23 { 24 Index pivot_row=j; 25 for(Index k=j+1;k<n;++k) 26 if(abs(A(k,j))>abs(A(pivot_row,j))) 27 pivot_row=j; 28 29 if(pivot_row!=j) 30 { 31 A.swap_rows(j,pivot_row); 32 std::swap(b(j),b(pivot_row)); 33 } 34 for(Index i=j+1;i<n;++i) 35 { 36 const double pivot=A(j,j); 37 if(pivot==0) 38 cerr<<"hehe"<<endl; 39 const double mult=A(i,j)/pivot; 40 A[i].slice(j)=scale_and_add(A[j].slice(j),-mult,A[i].slice(j)); 41 b(i)-=mult*b(j); 42 } 43 44 } 45 } 46 47 void classical_elimination(Matrix1&A,vector1&b) 48 { 49 const Index n=A.dim1(); 50 for(Index j=0;j<n-1;++j) 51 { 52 const double pivot=A(j,j); 53 if(pivot==0) 54 throw exception(); 55 56 for(Index i=j+1;i<n;++i) 57 { 58 const double mult =A(i,j)/pivot; 59 A[i].slice(j)=scale_and_add(A[j].slice(j),-mult,A[i].slice(j)); 60 b(i)-=mult*b(j); 61 } 62 } 63 } 64 65 vector1 back_substitution(const Matrix1&A,const vector1&b) 66 { 67 const Index n=A.dim1(); 68 vector1 x(n); 69 for(Index i=n-1;i>=0;--i) 70 { 71 double s=b(i)-dot_product(A[i].slice(i+1),x.slice(i+1)); 72 if(double m=A(i,i)) 73 x(i)=s/m; 74 else 75 throw exception(); 76 } 77 return x; 78 } 79 80 vector1 classical_gaussian_elimination(Matrix1 A,vector1 b) 81 { 82 //classical_elimination(A,b); 83 elim_with_partial_pivot(A,b); 84 return back_substitution(A,b); 85 } 86 87 vector1 operator*(const Matrix1&m,const vector1&u) 88 { 89 const Index n=m.dim1(); 90 vector1 v(n); 91 for(Index i=0;i<n;++i) 92 v(i)=dot_product(m[i],u); 93 return v; 94 } 95 96 void solve_random_system() 97 { 98 Matrix1 A(2,2); 99 A(0,0)=1; 100 A(0,1)=0; 101 A(1,0)=0; 102 A(1,1)=1; 103 vector1 b(2); 104 b(0)=5;b(1)=6; 105 106 cout<<A<<endl<<b<<endl; 107 try{ 108 vector1 x=classical_gaussian_elimination(A,b); 109 cout<<endl<<x<<endl; 110 vector1 v=A*x; 111 cout<<endl<<v<<endl; 112 } 113 catch(const exception&e) 114 { 115 cerr<<e.what()<<endl; 116 } 117 } 118 119 int main() 120 { 121 solve_random_system(); 122 123 while(1); 124 return 0; 125 126 }
习题4
1 template<class T,class P> 2 Matrix<T> apply1(P p,Matrix<T,1> a) 3 { 4 Matrix<T>t1=a; 5 for(Index i=0;i<t1.size();i++) 6 t1(i)=p(t1(i)); 7 return t1; 8 }
习题9
1 template<class T> 2 void swap_columns(Matrix<T,2>&m,Index i, Index j) 3 { 4 if (i == j) return; 5 6 T t; 7 int k=m.dim1(); 8 for(int h=0;h<k;h++) 9 { 10 t=m[h][i]; 11 m[h][i]=m[h][j]; 12 m[h][j]=t; 13 } 14 }
习题11
习题12 那个N我觉得得确定才能做= =求科普
时间: 2024-12-28 11:02:39