1! 2! CalculiX - A 3-dimensional finite element program 3! Copyright (C) 1998-2021 Guido Dhondt 4! 5! This program is free software; you can redistribute it and/or 6! modify it under the terms of the GNU General Public License as 7! published by the Free Software Foundation(version 2); 8! 9! 10! This program is distributed in the hope that it will be useful, 11! but WITHOUT ANY WARRANTY; without even the implied warranty of 12! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the 13! GNU General Public License for more details. 14! 15! You should have received a copy of the GNU General Public License 16! along with this program; if not, write to the Free Software 17! Foundation, Inc., 675 Mass Ave, Cambridge, MA 02139, USA. 18! 19 subroutine bounadd(node,is,ie,val,nodeboun,ndirboun,xboun, 20 & nboun,nboun_,iamboun,iamplitude,nam,ipompc,nodempc, 21 & coefmpc,nmpc,nmpc_,mpcfree,inotr,trab, 22 & ntrans,ikboun,ilboun,ikmpc,ilmpc,co,nk,nk_,labmpc,type, 23 & typeboun,nmethod,iperturb,fixed,vold,nodetrue,mi,label) 24! 25! adds a boundary condition to the data base 26! 27 implicit none 28! 29 logical fixed,rottracoupling 30! 31 character*1 type,typeboun(*) 32 character*20 labmpc(*),label 33! 34 integer nodeboun(*),ndirboun(*),node,is,ie,nboun,nboun_,i,j, 35 & iamboun(*),iamplitude,nam,ipompc(*),nodempc(3,*),nmpc,nmpc_, 36 & mpcfree,inotr(2,*),ntrans,ikboun(*),ilboun(*),ikmpc(*), 37 & ilmpc(*),itr,idof,newnode,number,id,idofnew,idnew,nk,nk_, 38 & mpcfreenew,nmethod,iperturb(*),ii,nodetrue,mi(*),three,kflag, 39 & iy(3),inumber,irotnode(11),irotdof(11) 40! 41 real*8 xboun(*),val,coefmpc(*),trab(7,*),a(3,3),co(3,*), 42 & vold(0:mi(2),*),dx(3) 43! 44 if(ntrans.le.0) then 45 itr=0 46 elseif(inotr(1,node).eq.0) then 47 itr=0 48 else 49 itr=inotr(1,node) 50 endif 51! 52! checking for boundary conditions on rotational dofs of 53! distributing couplings 54! 55 rottracoupling=.false. 56 if((ie.ge.4).and.(ie.le.6)) then 57! 58! rotational dof 59! 60 do ii=is,ie 61 irotnode(ii)=node 62 if(ii.gt.3) then 63c idof=8*(node-1)+ii+1 64 idof=8*(node-1)+ii 65 call nident(ikmpc,idof,nmpc,id) 66 if(id.gt.0) then 67 if(ikmpc(id).eq.idof) then 68 if(labmpc(ilmpc(id))(1:14).eq.'ROTTRACOUPLING')then 69 rottracoupling=.true. 70 irotnode(ii)= 71 & nodempc(1,nodempc(3,ipompc(ilmpc(id)))) 72 irotdof(ii)= 73 & nodempc(2,nodempc(3,ipompc(ilmpc(id)))) 74 itr=0 75 endif 76 endif 77 endif 78 endif 79 enddo 80 endif 81! 82 loop: do ii=is,ie 83! 84! change: transformations on rotations are taken into account 85! by the normal of the mean rotation MPC, not by expanding the 86! MPC in Carthesian coordinates 87! 88 if((itr.eq.0).or.(ii.eq.0).or.(ii.gt.3)) then 89! 90! no transformation applies: simple SPC 91! 92 if(rottracoupling) then 93 node=irotnode(ii) 94 if(ii.gt.3) then 95 i=irotdof(ii) 96 else 97 i=ii 98 endif 99 else 100c if(ii.le.3) then 101 if(ii.le.6) then 102 i=ii 103c elseif(ii.eq.4) then 104c i=5 105c elseif(ii.eq.5) then 106c i=6 107c elseif(ii.eq.6) then 108c i=7 109 elseif(ii.eq.8) then 110 i=4 111 elseif(ii.eq.11) then 112 i=0 113 else 114 write(*,*) '*ERROR in bounadd: unknown DOF: ', 115 & ii 116 call exit(201) 117 endif 118 endif 119! 120 if((fixed).and.(i.lt.5)) then 121 val=vold(i,nodetrue) 122 elseif(fixed) then 123 write(*,*) '*ERROR in bounadd: parameter FIXED cannot' 124 write(*,*) ' be used for rotations' 125 call exit(201) 126 endif 127 idof=8*(node-1)+i 128 call nident(ikboun,idof,nboun,id) 129 if(id.gt.0) then 130 if(ikboun(id).eq.idof) then 131 j=ilboun(id) 132 if(typeboun(j).ne.type) cycle loop 133 xboun(j)=val 134 if(nam.gt.0) iamboun(j)=iamplitude 135 cycle loop 136 endif 137 endif 138 nboun=nboun+1 139 if(nboun.gt.nboun_) then 140 write(*,*) '*ERROR in bounadd: increase nboun_' 141 call exit(201) 142 endif 143 if((nmethod.eq.4).and.(iperturb(1).le.1)) then 144 write(*,*) '*ERROR in bounadd: in a modal dynamic step' 145 write(*,*) ' new SPCs are not allowed' 146 call exit(201) 147 endif 148 nodeboun(nboun)=node 149 ndirboun(nboun)=i 150 xboun(nboun)=val 151 typeboun(nboun)=type 152 if(nam.gt.0) iamboun(nboun)=iamplitude 153! 154! updating ikboun and ilboun 155! 156 do j=nboun,id+2,-1 157 ikboun(j)=ikboun(j-1) 158 ilboun(j)=ilboun(j-1) 159 enddo 160 ikboun(id+1)=idof 161 ilboun(id+1)=nboun 162 else 163! 164! transformation applies: SPC is MPC in global carthesian 165! coordinates 166! 167 call transformatrix(trab(1,itr),co(1,node),a) 168c if(ii.le.3) then 169 if(ii.le.6) then 170 i=ii 171c elseif(ii.eq.4) then 172c i=5 173c elseif(ii.eq.5) then 174c i=6 175c elseif(ii.eq.6) then 176c i=7 177 elseif(ii.eq.8) then 178 i=4 179 elseif(ii.eq.11) then 180 i=0 181 else 182 write(*,*) '*ERROR in bounadd: unknown DOF: ', 183 & ii 184 call exit(201) 185 endif 186 if((fixed).and.(i.lt.5)) then 187 val=vold(i,nodetrue) 188 elseif(fixed) then 189 write(*,*) '*ERROR in bounadd: parameter FIXED cannot' 190 write(*,*) ' be used for rotations' 191 call exit(201) 192 endif 193 if(inotr(2,node).ne.0) then 194 newnode=inotr(2,node) 195 idofnew=8*(newnode-1)+i 196 call nident(ikboun,idofnew,nboun,idnew) 197 if(idnew.gt.0) then 198 if(ikboun(idnew).eq.idofnew) then 199 j=ilboun(idnew) 200c 201 if(typeboun(j).ne.type) cycle 202c 203 xboun(j)=val 204c typeboun(j)=type 205 if(nam.gt.0) iamboun(j)=iamplitude 206 cycle 207 endif 208 endif 209 else 210! 211! new node is generated for the inhomogeneous MPC term 212! 213 if((nmethod.eq.4).and.(iperturb(1).le.1)) then 214 write(*,*)'*ERROR in bounadd: in a modal dynamic step' 215 write(*,*) ' new SPCs are not allowed' 216 call exit(201) 217 endif 218 nk=nk+1 219 if(nk.gt.nk_) then 220 write(*,*) '*ERROR in bounadd: increase nk_' 221 call exit(201) 222 endif 223 newnode=nk 224 inotr(2,node)=newnode 225 idofnew=8*(newnode-1)+i 226 idnew=nboun 227! 228! copying the initial conditions from node into newnode 229! 230 do j=0,mi(2) 231 vold(j,newnode)=vold(j,node) 232 enddo 233 endif 234! 235! new mpc 236! 237 iy(1)=1 238 iy(2)=2 239 iy(3)=3 240 dx(1)=dabs(a(1,i)) 241 dx(2)=dabs(a(2,i)) 242 dx(3)=dabs(a(3,i)) 243 three=3 244 kflag=-2 245 call dsort(dx,iy,three,kflag) 246 do inumber=1,3 247 number=iy(inumber) 248 idof=8*(node-1)+number 249 call nident(ikmpc,idof,nmpc,id) 250 if(id.ne.0) then 251 if(ikmpc(id).eq.idof) cycle 252 endif 253 if(dabs(a(number,i)).lt.1.d-5) cycle 254 nmpc=nmpc+1 255 if(nmpc.gt.nmpc_) then 256 write(*,*) '*ERROR in bounadd: increase nmpc_' 257 call exit(201) 258 endif 259 labmpc(nmpc)=label 260 ipompc(nmpc)=mpcfree 261 do j=nmpc,id+2,-1 262 ikmpc(j)=ikmpc(j-1) 263 ilmpc(j)=ilmpc(j-1) 264 enddo 265 ikmpc(id+1)=idof 266 ilmpc(id+1)=nmpc 267 exit 268 enddo 269! 270! check whether a dependent term was found; if none was 271! found this can be due to the fact that: 272! - all dofs were used by other MPC's 273! - the MPC coefficients were too small 274! - or a combination of both 275! 276 if(inumber.gt.3) then 277 write(*,*) '*ERROR in bounadd' 278 write(*,*) ' SPC in node',node 279 write(*,*) ' and local direction',ii 280 write(*,*) ' cannot be applied: all' 281 write(*,*) ' degrees of freedom have' 282 write(*,*) ' been used by other MPCs' 283 write(*,*) ' or the coefficient is' 284 write(*,*) ' too small' 285 call exit(201) 286 endif 287! 288 inumber=inumber-1 289 do j=1,3 290 inumber=inumber+1 291 if(inumber.gt.3) inumber=1 292 number=iy(inumber) 293c if(dabs(a(number,i)).lt.1.d-5) cycle 294 if(dabs(a(number,i)).lt.1.d-30) cycle 295 nodempc(1,mpcfree)=node 296 nodempc(2,mpcfree)=number 297 coefmpc(mpcfree)=a(number,i) 298 mpcfree=nodempc(3,mpcfree) 299 if(mpcfree.eq.0) then 300 write(*,*) '*ERROR in bounadd: increase memmpc_' 301 call exit(201) 302 endif 303 enddo 304 nodempc(1,mpcfree)=newnode 305 nodempc(2,mpcfree)=i 306 coefmpc(mpcfree)=-1.d0 307 mpcfreenew=nodempc(3,mpcfree) 308 if(mpcfreenew.eq.0) then 309 write(*,*) '*ERROR in bounadd: increase nmpc_' 310 call exit(201) 311 endif 312 nodempc(3,mpcfree)=0 313 mpcfree=mpcfreenew 314! 315! nonhomogeneous term 316! 317 nboun=nboun+1 318 if(nboun.gt.nboun_) then 319 write(*,*) '*ERROR in bounadd: increase nboun_' 320 call exit(201) 321 endif 322 nodeboun(nboun)=newnode 323 ndirboun(nboun)=i 324 xboun(nboun)=val 325 typeboun(nboun)=type 326 if(nam.gt.0) iamboun(nboun)=iamplitude 327! 328! updating ikboun and ilboun 329! 330 do j=nboun,idnew+2,-1 331 ikboun(j)=ikboun(j-1) 332 ilboun(j)=ilboun(j-1) 333 enddo 334 ikboun(idnew+1)=idofnew 335 ilboun(idnew+1)=nboun 336! 337c enddo 338 endif 339 enddo loop 340! 341 return 342 end 343 344