Assuming you have the free dof numbers in d_free, the prescribed dofs in d_pres, the values of the non zero DBC in a_pres the solution to the free dofs a_free is: a_free = K[d_free, d_free] \ (f[d_free] - K[d_free, d_pres] * a_pres)
It will probably be faster than what you are doing now On Wednesday, August 24, 2016 at 6:29:52 AM UTC+2, Nguyen Vinh Phu wrote: > > Hello all, > > I am implementing a finite element solver in Julia. I have computed the > stiffness matrix (as a sparse matrix K) and the force vector (F). I have > some non zero boundary conditions on > the unknown. In Matlab, here is what I do: > > bcwt=mean(diag(K)); % a measure of the average size of an element in K > % used to keep the conditioning of the K matrix > udofs=fixedNode; % global indecies of the fixed displacements > f=f-K(:,udofs)*uFixed; % modify the force vector > K(udofs,:)=0; > K(:,udofs)=0; > K(udofs,udofs)=bcwt*speye(length(udofs)); % put ones*bcwt on the diagonal > f(udofs)=bcwt*speye(length(udofs))*uFixed; > U=K\f; > > I implemented the above in Julia and I got several errors associated with > the following > > f=f-K(:,udofs)*uFixed; % modify the force vector > K(udofs,:)=0; > K(:,udofs)=0; > K(udofs,udofs)=bcwt*speye(length(udofs)); % put ones*bcwt on the diagonal > f(udofs)=bcwt*speye(length(udofs))*uFixed; > > I would appreciate very much if someone can help. Thanks. > > Best regards, > Phu > >
