有限元刚度矩阵的一维变带宽存储用C++实现(一)

本文探讨了一维变带宽存储中辅助数组pDiag[n]的C++实现,用于存储有限元刚度矩阵的对角元地址。通过对节点自由度编号,计算单元定位矩阵,确定行最大半带宽,并提供程序流程图和代码示例。
摘要由CSDN通过智能技术生成
       有限元计算中的刚度矩阵是稀疏、对称矩阵,经过各个自由度的合理排序可以使得其中的非零元素集中在主对角线附近,形成带状矩阵。这样的稀疏带状矩阵在存储时一般有“二维等带宽”和“一维变带宽”存储两种方式。其中,“一维变带宽”存储可以节省大量的内存空间,提高计算速度,但现在网上大多数“一维变带宽”存储的资料都是基于 Fortran语言的,下面讨论“一维变带宽”存储的C++实现。以下仅仅是我自己学习的心得体会,如果大大们发现问题,欢迎提出讨论微笑
        限于篇幅,关于“一维变带宽”存储的详细原理就不赘述了,有兴趣的同学可以去看有限元程序设计的书籍,下面仅讨论其C++实现。

      “一维变带宽”存储需要一个存储刚度矩阵中带内元素的一维数组pGK[L](L为带内元素总数),和一个用于存储刚度矩阵中各行主对角线上元素在数组pGK[L]中的序号(称作对角元地址) 的一维辅助数组pDiag[n](n为刚度矩阵阶数)。由于刚度矩阵中元素关于主对角线对称,故pGK[L]中只存

有限元刚度矩阵一般都是一个稀疏矩阵,因此可以使用一维宽带储存来节省存储空间。下面是使用 MATLAB 求解一维宽带储存的有限元刚度矩阵的示例代码。 假设有一条长度为 L 的杆,将其分成 n 个小段,用线性元进行离散化,每个小段的长度为 h=L/n。杆的弹性模量为 E,截面积为 A,横向受力为 F。 首先定义杆的节点坐标和单元刚度矩阵: ``` n = 10; % 小段个数 h = L/n; % 小段长度 A = 1; % 杆的截面积 E = 1; % 杆的弹性模量 F = 1; % 杆的横向受力 % 定义节点坐标 X = linspace(0, L, n+1); % 定义单元刚度矩阵 k = [1 -1; -1 1] * A * E / h; ``` 然后根据单元刚度矩阵和节点坐标计算整个杆的刚度矩阵: ``` K = zeros(n+1, n+1); % 初始化刚度矩阵 for i = 1:n % 获取当前单元的节点编号 j = [i, i+1]; % 将当前单元的刚度矩阵加到整个杆的刚度矩阵上 K(j, j) = K(j, j) + k; end ``` 此时得到的刚度矩阵 K 是一个密集矩阵,需要进行压缩存储。使用一维宽带储存方法,可以将 K 压缩成一个一维数组: ``` % 获取非零元素的位置和值 [I, J, V] = find(K); nzmax = length(V); % 计算每一行的非零元素个数 nnz_row = sum(K~=0, 2); % 计算每个非零元素的偏移量 offset = zeros(n+1, 1); for i = 2:(n+1) offset(i) = offset(i-1) + nnz_row(i-1); end % 压缩存储刚度矩阵 K_sparse = zeros(1, nzmax + n + 1); K_sparse(1) = n+1; K_sparse(2:(n+2)) = offset'; K_sparse((n+3):end) = V'; ``` 最终得到的 K_sparse 就是一维宽带储存的刚度矩阵。其中 K_sparse(1) 表示矩阵的维数,K_sparse(2:n+2) 表示每一行非零元素的偏移量,K_sparse(n+3:end) 表示矩阵中的非零元素。
评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值