1 /* CalculiX - A 3-dimensional finite element program */
2 /* Copyright (C) 1998-2021 Guido Dhondt */
3
4 /* This program is free software; you can redistribute it and/or */
5 /* modify it under the terms of the GNU General Public License as */
6 /* published by the Free Software Foundation(version 2); */
7 /* */
8
9 /* This program is distributed in the hope that it will be useful, */
10 /* but WITHOUT ANY WARRANTY; without even the implied warranty of */
11 /* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the */
12 /* GNU General Public License for more details. */
13
14 /* You should have received a copy of the GNU General Public License */
15 /* along with this program; if not, write to the Free Software */
16 /* Foundation, Inc., 675 Mass Ave, Cambridge, MA 02139, USA. */
17
18 #include <stdlib.h>
19 #include <math.h>
20 #include <stdio.h>
21 #include <string.h>
22 #include "CalculiX.h"
23
24 #define min(a,b) ((a) <= (b) ? (a) : (b))
25 #define max(a,b) ((a) >= (b) ? (a) : (b))
26
mastructnmatrix(ITG * icols,ITG * jqs,ITG ** mast1p,ITG ** irowsp,ITG * ipointer,ITG * nzss,ITG * nactive,ITG * nnlconst)27 void mastructnmatrix(ITG *icols,ITG *jqs,ITG **mast1p,ITG **irowsp,
28 ITG *ipointer,ITG *nzss,ITG *nactive,ITG *nnlconst){
29
30 /* determines the structure of the (N^T*N)-Matrix;
31 (i.e. the location of the nonzeros */
32
33 ITG i,j,ii,jj,kk,index,jdof2,jdof1,nmast,ifree,kflag,indexe,isize,
34 *mast1=NULL,*irows=NULL,*next=NULL,jstart;
35
36 /* the indices in the comments follow FORTRAN convention, i.e. the
37 fields start with 1 */
38
39 mast1=*mast1p;
40 irows=*irowsp;
41 ifree=0;
42 kflag=2;
43
44 NNEW(next,ITG,*nzss);
45
46 for(kk=0;kk<*nnlconst;kk++){
47 jdof1=kk+1;
48
49 for(jj=0;jj<*nactive;jj++){
50 jdof2=jj+1;
51 if(jdof2>jdof1){
52 insert(ipointer,&mast1,&next,&jdof1,&jdof2,&ifree,nzss);
53 }
54 }
55 }
56
57 /* determination of the following fields:
58
59 - irows: row numbers, column per column
60 - jqs(i)= location in field irows of the first SUBdiagonal
61 nonzero in column i */
62
63 RENEW(irows,ITG,ifree);
64 nmast=0;
65 jqs[0]=1;
66 for(i=0;i<*nactive;i++){
67 index=ipointer[i];
68 do{
69 if(index==0) break;
70 irows[nmast++]=mast1[index-1];
71 index=next[index-1];
72 }while(1);
73 jqs[i+1]=nmast+1;
74 }
75
76 /* sorting the row numbers within each column */
77
78 for(i=0;i<*nactive;++i){
79 if(jqs[i+1]-jqs[i]>0){
80 isize=jqs[i+1]-jqs[i];
81 FORTRAN(isortii,(&irows[jqs[i]-1],&mast1[jqs[i]-1],&isize,&kflag));
82 }
83 }
84
85 /* removing duplicate entries */
86
87 nmast=0;
88 for(i=0;i<*nactive;i++){
89 jstart=nmast+1;
90 if(jqs[i+1]-jqs[i]>0){
91 irows[nmast++]=irows[jqs[i]-1];
92 for(j=jqs[i];j<jqs[i+1]-1;j++){
93 if(irows[j]==irows[nmast-1])continue;
94 irows[nmast++]=irows[j];
95 }
96 }
97 jqs[i]=jstart;
98 }
99 jqs[*nactive]=nmast+1;
100
101 for(i=0;i<*nactive;i++){
102 icols[i]=jqs[i+1]-jqs[i];
103 }
104 *nzss=jqs[*nactive]-1;
105
106 SFREE(next);
107
108 *mast1p=mast1;
109 *irowsp=irows;
110
111 return;
112
113 }
114