【刘庆源码共享】稀疏线性系统求解算法MGMRES(m) 之 矩阵类定义二(C++)

/* 
 * Copyright (c) 2009 湖南师范大学数计院 一心飞翔项目组
 * All Right Reserved
 *
 * 文件名：matrix.cpp  定义Point、Node、Matrix类的各个方法
 * 摘  要：定义矩阵类，包括矩阵的相关信息和方法
 *
 * 作  者：刘 庆
 * 修改日期：2009年7月19日21:15:12
 *
*/

 

/* 该类方法均基于稀疏线性系统的MCRF存储结构，MCRF相关介绍参见：http://blog.csdn.net/lq_0402/archive/2010/11/11/6003822.aspx */

 

/* 行初等变换 第base行乘以gene, 没有过滤掉乘以gene之后小于门槛值的元素 */
void Matrix::LnTransformMultiply(const long base,const Double& gene)
{
 long offset=data[base].index;
 /* 该行的非零非对角元素乘以gene */
 long i=0;
 for(i=data[base].index; i<data[base+1].index; i++)
 {
  data[i].value *= gene;
   data[offset].value = data[i].value;
   data[offset].index = data[i].index;
   offset++;
 }
 if(offset!=i)
 {
  for( ; i<total_elem+1; i++)
  {
   data[offset].value = data[i].value;
   data[offset].index = data[i].index;
   offset++;
  }
  total_elem = offset-1;
  for(long j=base+1; j<total_ln+1; j++)
   data[j].index -= i-offset;
 }

 /* 该行对角元素乘以gene */
 data[base].value *= gene;
}

/* 行初等变换 第base行加上第ln行  */
/* 算法：
 1、开辟一个辅助数组tempData存放一行的数据，其长度为total_col
 2、将base行的数据根据下标放于tempData中。
 3、将ln行的数据根据下标，加到tempData相应下标位置。
 4、开辟一个辅助数组tempDataAll，存放base和ln相加之后的整个新data的值
 5、将tempDataAll覆盖原有的data
*/
void Matrix::LnTransfromAdd(const long base,const long ln)
{
 Double *tempData = new Double[(total_col>total_ln?total_col:total_ln)+10];
 long nonZero = 0; // 非零非对角元素个素 */
 long i=0;
 for( i=0; i<total_col+1; i++)
  tempData[i] = 0;
 /* 将base行的元素置于tempData中 */
 /* 处理非零非对角元素 */
 for(i=data[base].index; i<data[base+1].index; i++ )
 {
  if(data[i].value!=0.0f)
  {
   tempData[data[i].index] = data[i].value;
   nonZero++;
  }
 }
 /* 处理对角元素 */
 tempData[base] = data[base].value;

 /* 将ln行的元素加到tempData相应位置 */
 /* 处理非零非对角元素 */
 for(i=data[ln].index; i<data[ln+1].index; i++)
 {
  // 先判断其值是否为0
  if(tempData[data[i].index]==0&&data[i].index!=base)
   nonZero++;
  tempData[data[i].index] += data[i].value;
  if(tempData[data[i].index]==0&&data[i].index!=base)
   nonZero--;
 }
 /* 处理对角元素 */
 if(tempData[ln]==0)
  nonZero++;
 tempData[ln] += data[ln].value;
 if(tempData[ln]==0)
  nonZero--;

 /* 将数组tempData中的数复制到data中的相应位置，涉及到三个步骤：
  1、新增一辅助data空间 tempDataAll
  2、将data和tempData中的数据复制到tempDataAll中
  3、用tempDataAll覆盖data
 */
 long overNum = nonZero-(data[base+1].index-data[base].index);
 total_elem = total_elem+overNum;
 Node tempDataAll;
 tempDataAll.GetSpace(total_elem+1);
 /* 将辅助数组初始化 */
 for(i=0; i<total_elem+1; i++)
 {
  tempDataAll[i].value = 0;
  tempDataAll[i].index = -1;
 }

 /* 前total_ln+1个位置 */
 for(i=0; i<total_ln+1; i++)
 {
  tempDataAll[i].value = data[i].value;
  if(i<=base)
   tempDataAll[i].index = data[i].index;
  else
  {
   tempDataAll[i].index = data[i].index+overNum;
  }
 }
 /* tota_ln+2 到 data[base].index-1 个位置 */
 for( ; i<data[base].index; i++){
  tempDataAll[i].value = data[i].value;
  tempDataAll[i].index = data[i].index;
 }
 long j=0;
 /* 新的 base 行的数据 */
 for( ; j<total_col; j++)
 {
  if(j==base){
   tempDataAll[base].value = tempData[j];
  }
  else if(tempData[j]!=0.0f)
  {
   tempDataAll[i].value = tempData[j];
   tempD