OpenCMISS-Iron Internal API Documentation
boundary_condition_routines.f90
Go to the documentation of this file.
1 
43 
45 MODULE boundary_conditions_routines
46 
47  USE base_routines
48  USE basis_routines
49  USE cmiss_mpi
50  USE comp_environment
51  USE constants
52  USE coordinate_routines
53  USE distributed_matrix_vector
54  USE domain_mappings
55  USE equations_set_constants
56  USE interface_conditions_constants
57  USE field_routines
58  USE input_output
59  USE iso_varying_string
60  USE kinds
61 #ifndef NOMPIMOD
62  USE mpi
63 #endif
64  USE node_routines
65  USE strings
66  USE timer
67  USE types
68  USE lists
69  USE linkedlist_routines
70 
71 #include "macros.h"
72 
73  IMPLICIT NONE
74 
75 #ifdef NOMPIMOD
76 #include "mpif.h"
77 #endif
78 
79  PRIVATE
80 
81  !Module parameters
82 
86  INTEGER(INTG), PARAMETER :: boundary_condition_dof_free=0
87  INTEGER(INTG), PARAMETER :: boundary_condition_dof_fixed=1
88  INTEGER(INTG), PARAMETER :: boundary_condition_dof_mixed=2
89  INTEGER(INTG), PARAMETER :: boundary_condition_dof_constrained=3
94  INTEGER(INTG), PARAMETER :: boundary_condition_free=0
95  INTEGER(INTG), PARAMETER :: boundary_condition_fixed=1
96  INTEGER(INTG), PARAMETER :: boundary_condition_fixed_inlet=2
97  INTEGER(INTG), PARAMETER :: boundary_condition_fixed_outlet=3
98  INTEGER(INTG), PARAMETER :: boundary_condition_fixed_wall=4
99  INTEGER(INTG), PARAMETER :: boundary_condition_moved_wall=5
100  INTEGER(INTG), PARAMETER :: boundary_condition_free_wall=6
101  INTEGER(INTG), PARAMETER :: boundary_condition_neumann_point=8
102  INTEGER(INTG), PARAMETER :: boundary_condition_neumann_integrated=9
103  INTEGER(INTG), PARAMETER :: boundary_condition_dirichlet=10
104  INTEGER(INTG), PARAMETER :: boundary_condition_cauchy=11
105  INTEGER(INTG), PARAMETER :: boundary_condition_robin=12
106  INTEGER(INTG), PARAMETER :: boundary_condition_fixed_incremented=13
107  INTEGER(INTG), PARAMETER :: boundary_condition_pressure=14
108  INTEGER(INTG), PARAMETER :: boundary_condition_pressure_incremented=15
109  INTEGER(INTG), PARAMETER :: boundary_condition_moved_wall_incremented=17
110  INTEGER(INTG), PARAMETER :: boundary_condition_correction_mass_increase=18
111  INTEGER(INTG), PARAMETER :: boundary_condition_impermeable_wall=19
112  INTEGER(INTG), PARAMETER :: boundary_condition_neumann_integrated_only=20
113  INTEGER(INTG), PARAMETER :: boundary_condition_linear_constraint=21
114  INTEGER(INTG), PARAMETER :: boundary_condition_neumann_point_incremented=22
115  INTEGER(INTG), PARAMETER :: boundary_condition_fixed_fitted=23
116  INTEGER(INTG), PARAMETER :: boundary_condition_fixed_nonreflecting=24
117  INTEGER(INTG), PARAMETER :: boundary_condition_fixed_cellml=25
118  INTEGER(INTG), PARAMETER :: boundary_condition_fixed_stree=26
120 
121  INTEGER(INTG), PARAMETER :: max_boundary_condition_number=26 !The maximum boundary condition type identifier, used for allocating an array with an entry for each type
122 
126  INTEGER(INTG), PARAMETER :: boundary_condition_sparse_matrices=1
127  INTEGER(INTG), PARAMETER :: boundary_condition_full_matrices=2
129 
130  !Module types
131 
132  !Module variables
133 
134  !Interfaces
135 
137  INTERFACE boundary_conditions_add_local_dof
138  MODULE PROCEDURE boundary_conditions_add_local_dof1
139  MODULE PROCEDURE boundary_conditions_add_local_dofs
140  END INTERFACE !BOUNDARY_CONDITIONS_ADD_LOCAL_DOF
141 
143  INTERFACE boundary_conditions_set_local_dof
144  MODULE PROCEDURE boundary_conditions_set_local_dof1
145  MODULE PROCEDURE boundary_conditions_set_local_dofs
146  END INTERFACE !BOUNDARY_CONDITIONS_SET_LOCAL_DOF
147 
148  PUBLIC boundary_condition_dof_free,boundary_condition_dof_fixed,boundary_condition_dof_mixed,boundary_condition_dof_constrained
149 
150  PUBLIC boundary_condition_free,boundary_condition_fixed,boundary_condition_fixed_inlet,&
151  & boundary_condition_fixed_outlet,boundary_condition_fixed_wall,boundary_condition_moved_wall,boundary_condition_free_wall,&
152  & boundary_condition_neumann_integrated,boundary_condition_dirichlet,boundary_condition_neumann_point, &
153  & boundary_condition_cauchy,boundary_condition_robin,boundary_condition_fixed_incremented,boundary_condition_pressure,&
154  & boundary_condition_pressure_incremented,boundary_condition_moved_wall_incremented, &
155  & boundary_condition_correction_mass_increase,boundary_condition_impermeable_wall,boundary_condition_neumann_integrated_only, &
156  & boundary_condition_neumann_point_incremented,boundary_condition_fixed_stree, &
157  & boundary_condition_fixed_fitted,boundary_condition_fixed_nonreflecting,boundary_condition_fixed_cellml
158 
159  PUBLIC boundary_condition_sparse_matrices,boundary_condition_full_matrices
160 
161  PUBLIC boundary_conditions_create_finish,boundary_conditions_create_start,boundary_conditions_destroy
162 
163  PUBLIC boundary_conditions_add_constant,boundary_conditions_add_local_dof,boundary_conditions_add_element, &
164  & boundary_conditions_add_node,boundary_conditions_variable_get
165 
166  PUBLIC boundary_conditions_set_constant,boundary_conditions_set_local_dof,boundary_conditions_set_element, &
167  & boundary_conditions_set_node,boundaryconditions_neumannintegrate,boundaryconditions_neumannsparsitytypeset
168 
169  PUBLIC boundaryconditions_constrainnodedofsequal
170 
171 CONTAINS
172 
173  !
174  !================================================================================================================================
175  !
176 
178  SUBROUTINE boundary_conditions_create_finish(BOUNDARY_CONDITIONS,ERR,ERROR,*)
179 
180  !Argument variables
181  TYPE(boundary_conditions_type), POINTER :: BOUNDARY_CONDITIONS
182  INTEGER(INTG), INTENT(OUT) :: ERR
183  TYPE(varying_string), INTENT(OUT) :: ERROR
184  !Local Variables
185  INTEGER(INTG) :: MPI_IERROR,SEND_COUNT,STORAGE_TYPE, NUMBER_OF_NON_ZEROS, NUMBER_OF_ROWS,COUNT
186  INTEGER(INTG) :: variable_idx,dof_idx, equ_matrix_idx, dirichlet_idx, row_idx, DUMMY, LAST, DIRICHLET_DOF
187  INTEGER(INTG) :: col_idx,equations_set_idx,parameterSetIdx
188  INTEGER(INTG) :: pressureIdx,neumannIdx
189  INTEGER(INTG), POINTER :: ROW_INDICES(:), COLUMN_INDICES(:)
190  TYPE(boundary_conditions_variable_type), POINTER :: BOUNDARY_CONDITION_VARIABLE
191  TYPE(domain_mapping_type), POINTER :: VARIABLE_DOMAIN_MAPPING
192  TYPE(field_variable_type), POINTER :: FIELD_VARIABLE
193  TYPE(boundary_conditions_dirichlet_type), POINTER :: BOUNDARY_CONDITIONS_DIRICHLET
194  TYPE(boundary_conditions_pressure_incremented_type), POINTER :: BOUNDARY_CONDITIONS_PRESSURE_INCREMENTED
195  TYPE(varying_string) :: LOCAL_ERROR
196  TYPE(solver_equations_type), POINTER :: SOLVER_EQUATIONS
197  TYPE(equations_set_type), POINTER :: EQUATIONS_SET
198  TYPE(equations_type), POINTER :: EQUATIONS
199  TYPE(equations_matrices_type), POINTER :: EQUATIONS_MATRICES
200  TYPE(equations_matrices_linear_type), POINTER :: LINEAR_MATRICES
201  TYPE(equations_matrices_dynamic_type), POINTER :: DYNAMIC_MATRICES
202  TYPE(equations_matrix_type), POINTER :: EQUATION_MATRIX
203  TYPE(boundary_conditions_sparsity_indices_type), POINTER :: SPARSITY_INDICES
204  TYPE(list_type), POINTER :: SPARSE_INDICES
205  TYPE(linkedlist),POINTER :: LIST(:)
206  INTEGER(INTG), ALLOCATABLE:: COLUMN_ARRAY(:)
207 
208  enters("BOUNDARY_CONDITIONS_CREATE_FINISH",err,error,*999)
209 
210  NULLIFY(boundary_conditions_pressure_incremented)
211 
212  IF(ASSOCIATED(boundary_conditions)) THEN
213  IF(boundary_conditions%BOUNDARY_CONDITIONS_FINISHED) THEN
214  CALL flagerror("Boundary conditions have already been finished.",err,error,*999)
215  ELSE
216  IF(ALLOCATED(boundary_conditions%BOUNDARY_CONDITIONS_VARIABLES)) THEN
217  IF(computational_environment%NUMBER_COMPUTATIONAL_NODES>0) THEN
218  !Transfer all the boundary conditions to all the computational nodes.
219  !\todo Look at this.
220  DO variable_idx=1,boundary_conditions%NUMBER_OF_BOUNDARY_CONDITIONS_VARIABLES
221  boundary_condition_variable=>boundary_conditions%BOUNDARY_CONDITIONS_VARIABLES(variable_idx)%PTR
222  IF(ASSOCIATED(boundary_condition_variable)) THEN
223  field_variable=>boundary_condition_variable%VARIABLE
224  IF(ASSOCIATED(field_variable)) THEN
225  variable_domain_mapping=>field_variable%DOMAIN_MAPPING
226  IF(ASSOCIATED(variable_domain_mapping)) THEN
227  send_count=variable_domain_mapping%NUMBER_OF_GLOBAL
228  !\todo This operation is a little expensive as we are doing an unnecessary sum across all the ranks in order to combin
229  !\todo the data from each rank into all ranks. We will see how this goes for now.
230  IF(computational_environment%NUMBER_COMPUTATIONAL_NODES>1) THEN
231  CALL mpi_allreduce(mpi_in_place,boundary_condition_variable%DOF_TYPES, &
232  & send_count,mpi_integer,mpi_sum,computational_environment%MPI_COMM,mpi_ierror)
233 
234  CALL mpi_error_check("MPI_ALLREDUCE",mpi_ierror,err,error,*999)
235  CALL mpi_allreduce(mpi_in_place,boundary_condition_variable%CONDITION_TYPES, &
236  & send_count,mpi_integer,mpi_sum,computational_environment%MPI_COMM,mpi_ierror)
237  CALL mpi_error_check("MPI_ALLREDUCE",mpi_ierror,err,error,*999)
238  ENDIF !mpi_in_place bug workaround - only do this when num comp nodes > 1
239 
240  ELSE
241  local_error="Field variable domain mapping is not associated for variable type "// &
242  & trim(number_to_vstring(variable_idx,"*",err,error))//"."
243  CALL flagerror(local_error,err,error,*999)
244  ENDIF
245 
246  IF(computational_environment%NUMBER_COMPUTATIONAL_NODES>1) THEN
247 
248  ! Update the total number of boundary condition types by summing across all nodes
249  CALL mpi_allreduce(mpi_in_place,boundary_condition_variable%DOF_COUNTS, &
250  & max_boundary_condition_number,mpi_integer,mpi_sum,computational_environment%MPI_COMM,mpi_ierror)
251  CALL mpi_error_check("MPI_ALLREDUCE",mpi_ierror,err,error,*999)
252  CALL mpi_allreduce(mpi_in_place,boundary_condition_variable%NUMBER_OF_DIRICHLET_CONDITIONS, &
253  & 1,mpi_integer,mpi_sum,computational_environment%MPI_COMM,mpi_ierror)
254  CALL mpi_error_check("MPI_ALLREDUCE",mpi_ierror,err,error,*999)
255  ENDif!mpi_in_place bug workaround - only do this when num comp nodes > 1
256  ! Check that the boundary conditions set are appropriate for equations sets
257  CALL boundaryconditions_checkequations(boundary_condition_variable,err,error,*999)
258 
259  IF(computational_environment%NUMBER_COMPUTATIONAL_NODES>1) THEN
260  !Make sure the required parameter sets are created on all computational nodes and begin updating them
261  CALL mpi_allreduce(mpi_in_place,boundary_condition_variable%parameterSetRequired, &
262  & field_number_of_set_types,mpi_logical,mpi_lor,computational_environment%MPI_COMM,mpi_ierror)
263  CALL mpi_error_check("MPI_ALLREDUCE",mpi_ierror,err,error,*999)
264  DO parametersetidx=1,field_number_of_set_types
265  IF(boundary_condition_variable%parameterSetRequired(parametersetidx)) THEN
266  CALL field_parametersetensurecreated(field_variable%FIELD,field_variable%VARIABLE_TYPE, &
267  & parametersetidx,err,error,*999)
268  CALL field_parameter_set_update_start(field_variable%FIELD,field_variable%VARIABLE_TYPE, &
269  & parametersetidx,err,error,*999)
270  END IF
271  END DO
272  ENDif!mpi_in_place bug workaround - only do this when num comp nodes > 1
273  ! Set up pressure incremented condition, if it exists
274  IF(boundary_condition_variable%DOF_COUNTS(boundary_condition_pressure_incremented)>0) THEN
275  CALL boundary_conditions_pressure_incremented_initialise(boundary_condition_variable,err,error,*999)
276  boundary_conditions_pressure_incremented=>boundary_condition_variable%PRESSURE_INCREMENTED_BOUNDARY_CONDITIONS
277  END IF
278 
279  ! Set up Neumann condition information if there are any Neumann conditions
280  IF(boundary_condition_variable%DOF_COUNTS(boundary_condition_neumann_point)>0.OR. &
281  & boundary_condition_variable%DOF_COUNTS(boundary_condition_neumann_point_incremented)>0) THEN
282  CALL boundaryconditions_neumanninitialise(boundary_condition_variable,err,error,*999)
283  END IF
284 
285  ! Loop over all global DOFs, keeping track of the dof indices of specific BC types where required
286  pressureidx=1
287  neumannidx=1
288  DO dof_idx=1,field_variable%NUMBER_OF_GLOBAL_DOFS
289  IF(boundary_condition_variable%CONDITION_TYPES(dof_idx)== boundary_condition_pressure_incremented) THEN
290  boundary_conditions_pressure_incremented%PRESSURE_INCREMENTED_DOF_INDICES(pressureidx)=dof_idx
291  pressureidx=pressureidx+1
292  ELSE IF(boundary_condition_variable%CONDITION_TYPES(dof_idx)==boundary_condition_neumann_point.OR. &
293  & boundary_condition_variable%CONDITION_TYPES(dof_idx)==boundary_condition_neumann_point_incremented) THEN
294  boundary_condition_variable%neumannBoundaryConditions%setDofs(neumannidx)=dof_idx
295  neumannidx=neumannidx+1
296  END IF
297  END DO
298 
299  ! Now that we know where Neumann point DOFs are, we can calculate matrix structure
300  IF(boundary_condition_variable%DOF_COUNTS(boundary_condition_neumann_point)>0.OR. &
301  & boundary_condition_variable%DOF_COUNTS(boundary_condition_neumann_point_incremented)>0) THEN
302  CALL boundaryconditions_neumannmatricesinitialise(boundary_condition_variable,err,error,*999)
303  END IF
304 
305  ! Check that there is at least one dirichlet condition
306  IF(boundary_condition_variable%NUMBER_OF_DIRICHLET_CONDITIONS>0) THEN
307  CALL boundary_conditions_dirichlet_initialise(boundary_condition_variable,err,error,*999)
308  boundary_conditions_dirichlet=>boundary_condition_variable%DIRICHLET_BOUNDARY_CONDITIONS
309  IF(ASSOCIATED(boundary_conditions_dirichlet)) THEN
310  ! Find dirichlet conditions
311  dirichlet_idx=1
312  DO dof_idx=1,field_variable%NUMBER_OF_GLOBAL_DOFS
313  IF(boundary_condition_variable%DOF_TYPES(dof_idx)==boundary_condition_dof_fixed) THEN
314  boundary_conditions_dirichlet%DIRICHLET_DOF_INDICES(dirichlet_idx)=dof_idx
315  dirichlet_idx=dirichlet_idx+1
316  ENDIF
317  ENDDO
318 
319  !Store Dirichlet dof indices
320  solver_equations=>boundary_conditions%SOLVER_EQUATIONS
321  IF(ASSOCIATED(solver_equations)) THEN
322  IF(ASSOCIATED(solver_equations%SOLVER_MAPPING)) THEN
323  DO equations_set_idx=1,solver_equations%SOLVER_MAPPING%NUMBER_OF_EQUATIONS_SETS
324  equations_set=>solver_equations%SOLVER_MAPPING%EQUATIONS_SETS(equations_set_idx)%PTR
325  IF(ASSOCIATED(equations_set)) THEN
326  equations=>equations_set%EQUATIONS
327  IF(ASSOCIATED(equations)) THEN
328  equations_matrices=>equations%EQUATIONS_MATRICES
329  IF(ASSOCIATED(equations_matrices)) THEN
330  linear_matrices=>equations_matrices%LINEAR_MATRICES
331  IF(ASSOCIATED(linear_matrices)) THEN
332  !Iterate through equations matrices
333  DO equ_matrix_idx=1,linear_matrices%NUMBER_OF_LINEAR_MATRICES
334  equation_matrix=>linear_matrices%MATRICES(equ_matrix_idx)%PTR
335  CALL distributed_matrix_storage_type_get(equation_matrix%MATRIX,storage_type,err,error,*999)
336  IF(ASSOCIATED(equation_matrix)) THEN
337  SELECT CASE(storage_type)
338  CASE(distributed_matrix_block_storage_type)
339  !Do nothing
340  CASE(distributed_matrix_diagonal_storage_type)
341  !Do nothing
342  CASE(distributed_matrix_column_major_storage_type)
343  CALL flagerror("Not implemented for column major storage.",err,error,*999)
344  CASE(distributed_matrix_row_major_storage_type)
345  CALL flagerror("Not implemented for row major storage.",err,error,*999)
346  CASE(distributed_matrix_compressed_row_storage_type)
347  !Get Sparsity pattern, number of non zeros, number of rows
348  CALL distributed_matrix_storage_locations_get(equation_matrix%MATRIX,row_indices, &
349  & column_indices,err,error,*999)
350  CALL distributed_matrix_number_non_zeros_get(equation_matrix%MATRIX,number_of_non_zeros, &
351  & err,error,*999)
352  !Get the matrix stored as a linked list
353  CALL distributed_matrix_linklist_get(equation_matrix%MATRIX,list,err,error,*999)
354  number_of_rows=equations_matrices%TOTAL_NUMBER_OF_ROWS
355  !Initialise sparsity indices arrays
356  CALL boundaryconditions_sparsityindicesinitialise(boundary_conditions_dirichlet% &
357  & linear_sparsity_indices(equations_set_idx,equ_matrix_idx)%PTR, &
358  & boundary_condition_variable%NUMBER_OF_DIRICHLET_CONDITIONS,err,error,*999)
359  !Find dirichlet columns and store the non zero indices (with respect to the 1D storage array)
360  NULLIFY(sparsity_indices)
361  sparsity_indices=>boundary_conditions_dirichlet%LINEAR_SPARSITY_INDICES( &
362  & equations_set_idx,equ_matrix_idx)%PTR
363  IF(ASSOCIATED(sparsity_indices)) THEN
364  !Setup list for storing dirichlet non zero indices
365  NULLIFY(sparse_indices)
366  CALL list_create_start(sparse_indices,err,error,*999)
367  CALL list_data_type_set(sparse_indices,list_intg_type,err,error,*999)
368  CALL list_initial_size_set(sparse_indices, &
369  & boundary_condition_variable%NUMBER_OF_DIRICHLET_CONDITIONS*( &
370  & number_of_non_zeros/number_of_rows),err,error,*999)
371  CALL list_create_finish(sparse_indices,err,error,*999)
372  count=0
373  sparsity_indices%SPARSE_COLUMN_INDICES(1)=1
374  last=1
375  DO dirichlet_idx=1,boundary_condition_variable%NUMBER_OF_DIRICHLET_CONDITIONS
376  dirichlet_dof=boundary_conditions_dirichlet%DIRICHLET_DOF_INDICES(dirichlet_idx)
377  CALL linkedlist_to_array(list(dirichlet_dof),column_array,err,error,*999)
378  DO row_idx=1,size(column_array)
379  CALL list_item_add(sparse_indices,column_array(row_idx),err,error,*999)
380  count=count+1
381  last=row_idx+1
382  ENDDO
383  sparsity_indices%SPARSE_COLUMN_INDICES(dirichlet_idx+1)=count+1
384  ENDDO
385  CALL list_detach_and_destroy(sparse_indices,dummy,sparsity_indices%SPARSE_ROW_INDICES, &
386  & err,error,*999)
387  DO col_idx =1,number_of_rows
388  CALL linkedlist_destroy(list(col_idx),err,error,*999)
389  ENDDO
390  ELSE
391  local_error="Sparsity indices arrays are not associated for this equations matrix."
392  CALL flagerror(local_error,err,error,*999)
393  ENDIF
394  CASE(distributed_matrix_compressed_column_storage_type)
395  CALL flagerror("Not implemented for compressed column storage.",err,error,*999)
396  CASE(distributed_matrix_row_column_storage_type)
397  CALL flagerror("Not implemented for row column storage.",err,error,*999)
398  CASE DEFAULT
399  local_error="The storage type of "//trim(number_to_vstring(storage_type,"*",err,error)) &
400  //" is invalid."
401  CALL flagerror(local_error,err,error,*999)
402  END SELECT
403  ELSE
404  CALL flagerror("The equation matrix is not associated.",err,error,*999)
405  ENDIF
406  ENDDO
407  ENDIF
408 
409  dynamic_matrices=>equations_matrices%DYNAMIC_MATRICES
410  IF(ASSOCIATED(dynamic_matrices)) THEN
411  !Iterate through equations matrices
412  DO equ_matrix_idx=1,dynamic_matrices%NUMBER_OF_DYNAMIC_MATRICES
413  equation_matrix=>dynamic_matrices%MATRICES(equ_matrix_idx)%PTR
414  CALL distributed_matrix_storage_type_get(equation_matrix%MATRIX,storage_type,err,error,*999)
415  IF(ASSOCIATED(equation_matrix)) THEN
416  SELECT CASE(storage_type)
417  CASE(distributed_matrix_block_storage_type)
418  !Do nothing
419  CASE(distributed_matrix_diagonal_storage_type)
420  !Do nothing
421  CASE(distributed_matrix_column_major_storage_type)
422  CALL flagerror("Not implemented for column major storage.",err,error,*999)
423  CASE(distributed_matrix_row_major_storage_type)
424  CALL flagerror("Not implemented for row major storage.",err,error,*999)
425  CASE(distributed_matrix_compressed_row_storage_type)
426  !Get Sparsity pattern, number of non zeros, number of rows
427  CALL distributed_matrix_storage_locations_get(equation_matrix%MATRIX,row_indices, &
428  & column_indices,err,error,*999)
429  CALL distributed_matrix_number_non_zeros_get(equation_matrix%MATRIX,number_of_non_zeros, &
430  & err,error,*999)
431  !Sparse matrix in a list
432  CALL distributed_matrix_linklist_get(equation_matrix%MATRIX,list,err,error,*999)
433  number_of_rows=equations_matrices%TOTAL_NUMBER_OF_ROWS
434  !Intialise sparsity indices arrays
435  CALL boundaryconditions_sparsityindicesinitialise(boundary_conditions_dirichlet% &
436  & dynamic_sparsity_indices(equations_set_idx,equ_matrix_idx)%PTR, &
437  & boundary_condition_variable%NUMBER_OF_DIRICHLET_CONDITIONS,err,error,*999)
438  !Find dirichlet columns and store the non zero indices (with respect to the 1D storage array)
439  NULLIFY(sparsity_indices)
440  sparsity_indices=>boundary_conditions_dirichlet%DYNAMIC_SPARSITY_INDICES( &
441  & equations_set_idx,equ_matrix_idx)%PTR
442  IF(ASSOCIATED(sparsity_indices)) THEN
443  ! Setup list for storing dirichlet non zero indices
444  NULLIFY(sparse_indices)
445  CALL list_create_start(sparse_indices,err,error,*999)
446  CALL list_data_type_set(sparse_indices,list_intg_type,err,error,*999)
447  CALL list_initial_size_set(sparse_indices, &
448  & boundary_condition_variable%NUMBER_OF_DIRICHLET_CONDITIONS*( &
449  & number_of_non_zeros/number_of_rows),err,error,*999)
450  CALL list_create_finish(sparse_indices,err,error,*999)
451  count=0
452  sparsity_indices%SPARSE_COLUMN_INDICES(1)=1
453  last=1
454  DO dirichlet_idx=1,boundary_condition_variable%NUMBER_OF_DIRICHLET_CONDITIONS
455  !Dirichlet columns
456  dirichlet_dof=boundary_conditions_dirichlet%DIRICHLET_DOF_INDICES(dirichlet_idx)
457  CALL linkedlist_to_array(list(dirichlet_dof),column_array,err,error,*999)
458  !The row indices
459  DO row_idx=1,size(column_array)
460  CALL list_item_add(sparse_indices,column_array(row_idx),err,error,*999)
461  count=count+1
462  last=row_idx+1
463  ENDDO
464  sparsity_indices%SPARSE_COLUMN_INDICES(dirichlet_idx+1)=count+1
465  ENDDO
466  CALL list_detach_and_destroy(sparse_indices,dummy,sparsity_indices%SPARSE_ROW_INDICES, &
467  & err,error,*999)
468  DO col_idx =1,number_of_rows
469  CALL linkedlist_destroy(list(col_idx),err,error,*999)
470  ENDDO
471  ELSE
472  local_error="Sparsity indices arrays are not associated for this equations matrix."
473  CALL flagerror(local_error,err,error,*999)
474  ENDIF
475  CASE(distributed_matrix_compressed_column_storage_type)
476  CALL flagerror("Not implemented for compressed column storage.",err,error,*999)
477  CASE(distributed_matrix_row_column_storage_type)
478  CALL flagerror("Not implemented for row column storage.",err,error,*999)
479  CASE DEFAULT
480  local_error="The storage type of "//trim(number_to_vstring(storage_type,"*",err,error)) &
481  //" is invalid."
482  CALL flagerror(local_error,err,error,*999)
483  END SELECT
484  ELSE
485  CALL flagerror("The equation matrix is not associated.",err,error,*999)
486  ENDIF
487  ENDDO
488  ENDIF
489  ELSE
490  local_error="Equations Matrices is not associated for these Equations."
491  CALL flagerror(local_error,err,error,*999)
492  ENDIF
493  ELSE
494  local_error="Equations is not associated for this Equations Set."
495  CALL flagerror(local_error,err,error,*999)
496  ENDIF
497  ELSE
498  local_error="Equations Set is not associated for boundary conditions variable "// &
499  & trim(number_to_vstring(variable_idx,"*",err,error))//"."
500  CALL flagerror(local_error,err,error,*999)
501  ENDIF
502  ENDDO !equations_set_idx
503  !\todo Update interface sparsity structure calculate first then update code below.
504 ! !Loop over interface conditions. Note that only linear interface matrices implemented so far.
505 ! DO interface_condition_idx=1,SOLVER_EQUATIONS%SOLVER_MAPPING%NUMBER_OF_INTERFACE_CONDITIONS
506 ! INTERFACE_CONDITION=>SOLVER_EQUATIONS%SOLVER_MAPPING%INTERFACE_CONDITIONS(interface_condition_idx)%PTR
507 ! IF(ASSOCIATED(INTERFACE_CONDITION)) THEN
508 ! INTERFACE_EQUATIONS=>INTERFACE_CONDITION%INTERFACE_EQUATIONS
509 ! IF(ASSOCIATED(INTERFACE_EQUATIONS)) THEN
510 ! INTERFACE_MATRICES=>INTERFACE_EQUATIONS%INTERFACE_MATRICES
511 ! IF(ASSOCIATED(INTERFACE_MATRICES)) THEN
512 ! !Iterate through equations matrices
513 ! DO interface_matrix_idx=1,INTERFACE_MATRICES%NUMBER_OF_INTERFACE_MATRICES
514 ! INTERFACE_MATRIX=>INTERFACE_MATRICES%MATRICES(interface_matrix_idx)%PTR
515 ! IF(ASSOCIATED(INTERFACE_MATRIX)) THEN
516 ! CALL DISTRIBUTED_MATRIX_STORAGE_TYPE_GET(INTERFACE_MATRIX%MATRIX,STORAGE_TYPE,ERR,ERROR,*999)
517 ! SELECT CASE(STORAGE_TYPE)
518 ! CASE(DISTRIBUTED_MATRIX_BLOCK_STORAGE_TYPE)
519 ! !Do nothing
520 ! CASE(DISTRIBUTED_MATRIX_DIAGONAL_STORAGE_TYPE)
521 ! !Do nothing
522 ! CASE(DISTRIBUTED_MATRIX_COLUMN_MAJOR_STORAGE_TYPE)
523 ! CALL FlagError("Not implemented for column major storage.",ERR,ERROR,*999)
524 ! CASE(DISTRIBUTED_MATRIX_ROW_MAJOR_STORAGE_TYPE)
525 ! CALL FlagError("Not implemented for row major storage.",ERR,ERROR,*999)
526 ! CASE(DISTRIBUTED_MATRIX_COMPRESSED_ROW_STORAGE_TYPE)
527 ! !Get Sparsity pattern, number of non zeros, number of rows
528 ! CALL DISTRIBUTED_MATRIX_STORAGE_LOCATIONS_GET(INTERFACE_MATRIX%MATRIX,ROW_INDICES, &
529 ! & COLUMN_INDICES,ERR,ERROR,*999)
530 ! CALL DISTRIBUTED_MATRIX_NUMBER_NON_ZEROS_GET(INTERFACE_MATRIX%MATRIX,NUMBER_OF_NON_ZEROS, &
531 ! & ERR,ERROR,*999)
532 ! !Get the matrix stored as a linked list
533 ! CALL DISTRIBUTED_MATRIX_LINKLIST_GET(INTERFACE_MATRIX%MATRIX,LIST,ERR,ERROR,*999)
534 ! NUMBER_OF_ROWS=EQUATIONS_MATRICES%TOTAL_NUMBER_OF_ROWS
535 ! !Initialise sparsity indices arrays
536 ! CALL BoundaryConditions_SparsityIndicesInitialise(BOUNDARY_CONDITIONS_DIRICHLET% &
537 ! & LINEAR_SPARSITY_INDICES(interface_condition_idx,interface_matrix_idx)%PTR, &
538 ! & BOUNDARY_CONDITION_VARIABLE%NUMBER_OF_DIRICHLET_CONDITIONS,ERR,ERROR,*999)
539 ! !Find dirichlet columns and store the non zero indices (with respect to the 1D storage array)
540 ! NULLIFY(SPARSITY_INDICES)
541 ! SPARSITY_INDICES=>BOUNDARY_CONDITIONS_DIRICHLET%LINEAR_SPARSITY_INDICES( &
542 ! & interface_condition_idx,interface_matrix_idx)%PTR
543 ! IF(ASSOCIATED(SPARSITY_INDICES)) THEN
544 ! !Setup list for storing dirichlet non zero indices
545 ! NULLIFY(SPARSE_INDICES)
546 ! CALL LIST_CREATE_START(SPARSE_INDICES,ERR,ERROR,*999)
547 ! CALL LIST_DATA_TYPE_SET(SPARSE_INDICES,LIST_INTG_TYPE,ERR,ERROR,*999)
548 ! CALL LIST_INITIAL_SIZE_SET(SPARSE_INDICES, &
549 ! & NUMBER_OF_DIRICHLET_CONDITIONS*(NUMBER_OF_NON_ZEROS/NUMBER_OF_ROWS),ERR,ERROR,*999)
550 ! CALL LIST_CREATE_FINISH(SPARSE_INDICES,ERR,ERROR,*999)
551 ! COUNT=0
552 ! SPARSITY_INDICES%SPARSE_COLUMN_INDICES(1)=1
553 ! LAST=1
554 ! DO dirichlet_idx=1,BOUNDARY_CONDITION_VARIABLE%NUMBER_OF_DIRICHLET_CONDITIONS
555 ! DIRICHLET_DOF=BOUNDARY_CONDITIONS_DIRICHLET%DIRICHLET_DOF_INDICES(dirichlet_idx)
556 ! CALL LinkedList_to_Array(list(DIRICHLET_DOF),column_array)
557 ! DO row_idx=1,size(column_array)
558 ! CALL LIST_ITEM_ADD(SPARSE_INDICES,column_array(row_idx),ERR,ERROR,*999)
559 ! COUNT=COUNT+1
560 ! LAST=row_idx+1
561 ! ENDDO
562 ! SPARSITY_INDICES%SPARSE_COLUMN_INDICES(dirichlet_idx+1)=COUNT+1
563 ! ENDDO
564 ! CALL LIST_DETACH_AND_DESTROY(SPARSE_INDICES,DUMMY,SPARSITY_INDICES%SPARSE_ROW_INDICES, &
565 ! & ERR,ERROR,*999)
566 ! DO col_idx =1,NUMBER_OF_ROWS
567 ! CALL LINKEDLIST_DESTROY(list(col_idx))
568 ! ENDDO
569 ! ELSE
570 ! LOCAL_ERROR="Sparsity indices arrays are not associated for this interface matrix."
571 ! CALL FlagError(LOCAL_ERROR,ERR,ERROR,*999)
572 ! ENDIF
573 ! CASE(DISTRIBUTED_MATRIX_COMPRESSED_COLUMN_STORAGE_TYPE)
574 ! CALL FlagError("Not implemented for compressed column storage.",ERR,ERROR,*999)
575 ! CASE(DISTRIBUTED_MATRIX_ROW_COLUMN_STORAGE_TYPE)
576 ! CALL FlagError("Not implemented for row column storage.",ERR,ERROR,*999)
577 ! CASE DEFAULT
578 ! LOCAL_ERROR="The storage type of "//TRIM(NUMBER_TO_VSTRING(STORAGE_TYPE,"*",ERR,ERROR)) &
579 ! //" is invalid."
580 ! CALL FlagError(LOCAL_ERROR,ERR,ERROR,*999)
581 ! END SELECT
582 ! ELSE
583 ! CALL FlagError("The interface matrix is not associated.",ERR,ERROR,*999)
584 ! ENDIF
585 ! ENDDO
586 ! ELSE
587 ! LOCAL_ERROR="Interface matrices is not associated for these interface equations."
588 ! CALL FlagError(LOCAL_ERROR,ERR,ERROR,*999)
589 ! ENDIF
590 ! ELSE
591 ! LOCAL_ERROR="Interface equations is not associated for this interface condition."
592 ! CALL FlagError(LOCAL_ERROR,ERR,ERROR,*999)
593 ! ENDIF
594 ! ELSE
595 ! LOCAL_ERROR="Interface condition is not associated for boundary conditions variable "// &
596 ! & TRIM(NUMBER_TO_VSTRING(variable_idx,"*",ERR,ERROR))//"."
597 ! CALL FlagError(LOCAL_ERROR,ERR,ERROR,*999)
598 ! ENDIF
599 ! ENDDO !interface_condition_idx
600  ELSE
601  local_error="Solver equations solver mapping is not associated."
602  CALL flagerror(local_error,err,error,*999)
603  ENDIF
604  ELSE
605  local_error="Solver equations is not associated."
606  CALL flagerror(local_error,err,error,*999)
607  ENDIF
608  ELSE
609  local_error="Dirichlet Boundary Conditions type is not associated for boundary condition variable type "// &
610  & trim(number_to_vstring(variable_idx,"*",err,error))//"."
611  CALL flagerror(local_error,err,error,*999)
612  ENDIF
613  ENDIF
614  ! Finish field update
615  DO parametersetidx=1,field_number_of_set_types
616  IF(boundary_condition_variable%parameterSetRequired(parametersetidx)) THEN
617  CALL field_parameter_set_update_finish(field_variable%FIELD,field_variable%VARIABLE_TYPE, &
618  & parametersetidx,err,error,*999)
619  END IF
620  END DO
621 
622  !Finish creating the boundary conditions DOF constraints
623  CALL boundaryconditions_dofconstraintscreatefinish(boundary_condition_variable,err,error,*999)
624  ELSE
625  local_error="Field variable is not associated for variable index "// &
626  & trim(number_to_vstring(variable_idx,"*",err,error))//"."
627  CALL flagerror(local_error,err,error,*999)
628  ENDIF
629  ELSE
630  CALL flagerror("Boundary conditions variable is not associated for variable index "// &
631  & trim(number_to_vstring(variable_idx,"*",err,error))//".",err,error,*999)
632  ENDIF
633  ENDDO ! variable_idx
634 
635  ENDIF
636  !Set the finished flag
637  boundary_conditions%BOUNDARY_CONDITIONS_FINISHED=.true.
638  ELSE
639  CALL flagerror("Boundary conditions variables array is not allocated.",err,error,*999)
640  ENDIF
641  ENDIF
642  ELSE
643  CALL flagerror("Boundary conditions is not associated.",err,error,*999)
644  ENDIF
645  IF(diagnostics1) THEN
646  CALL write_string(diagnostic_output_type,"Boundary conditions:",err,error,*999)
647  DO variable_idx=1,boundary_conditions%NUMBER_OF_BOUNDARY_CONDITIONS_VARIABLES
648  boundary_condition_variable=>boundary_conditions%BOUNDARY_CONDITIONS_VARIABLES(variable_idx)%PTR
649  CALL write_string_value(diagnostic_output_type," Variable type = ",boundary_condition_variable%VARIABLE_TYPE, &
650  & err,error,*999)
651  IF(ASSOCIATED(boundary_condition_variable)) THEN
652  field_variable=>boundary_condition_variable%VARIABLE
653  variable_domain_mapping=>field_variable%DOMAIN_MAPPING
654  CALL write_string_value(diagnostic_output_type," Number of global dofs = ",variable_domain_mapping% &
655  & number_of_global,err,error,*999)
656  CALL write_string_vector(diagnostic_output_type,1,1,variable_domain_mapping%NUMBER_OF_GLOBAL,8,8, &
657  & boundary_condition_variable%CONDITION_TYPES,'(" Global BCs:",8(X,I8))','(15X,8(X,I8))', &
658  & err,error,*999)
659  ELSE
660  CALL flagerror("Boundary condition variable is not associated",err,error,*999)
661  ENDIF
662  ENDDO !variable_idx
663  ENDIF
664 
665  exits("BOUNDARY_CONDITIONS_CREATE_FINISH")
666  RETURN
667 999 errorsexits("BOUNDARY_CONDITIONS_CREATE_FINISH",err,error)
668  RETURN 1
669 
670  END SUBROUTINE boundary_conditions_create_finish
671 
672  !
673  !================================================================================================================================
674  !
675 
677  SUBROUTINE boundary_conditions_create_start(SOLVER_EQUATIONS,BOUNDARY_CONDITIONS,ERR,ERROR,*)
678 
679  !Argument variables
680  TYPE(solver_equations_type), POINTER :: SOLVER_EQUATIONS
681  TYPE(boundary_conditions_type), POINTER :: BOUNDARY_CONDITIONS
682  INTEGER(INTG), INTENT(OUT) :: ERR
683  TYPE(varying_string), INTENT(OUT) :: ERROR
684  !Local Variables
685  TYPE(varying_string) :: LOCAL_ERROR
686 
687  enters("BOUNDARY_CONDITIONS_CREATE_START",err,error,*999)
688 
689  IF(ASSOCIATED(solver_equations)) THEN
690  IF(ASSOCIATED(solver_equations%BOUNDARY_CONDITIONS)) THEN
691  CALL flagerror("Boundary conditions are already associated for the solver equations.",err,error,*999)
692  ELSE
693  IF(ASSOCIATED(boundary_conditions)) THEN
694  CALL flagerror("Boundary conditions is already associated.",err,error,*999)
695  ELSE
696  IF(ASSOCIATED(solver_equations%SOLVER_MAPPING)) THEN
697  !Initialise the boundary conditions
698  CALL boundary_conditions_initialise(solver_equations,err,error,*999)
699  ELSE
700  local_error="Solver equations solver mapping is not associated."
701  CALL flagerror(local_error,err,error,*999)
702  ENDIF
703  !Return the pointer
704  boundary_conditions=>solver_equations%BOUNDARY_CONDITIONS
705  ENDIF
706  ENDIF
707  ELSE
708  CALL flagerror("Solver equations is not associated.",err,error,*999)
709  ENDIF
710 
711  exits("BOUNDARY_CONDITIONS_CREATE_START")
712  RETURN
713 999 errorsexits("BOUNDARY_CONDITIONS_CREATE_START",err,error)
714  RETURN 1
715 
716  END SUBROUTINE boundary_conditions_create_start
717 
718  !
719  !================================================================================================================================
720  !
721 
723  SUBROUTINE boundary_conditions_destroy(BOUNDARY_CONDITIONS,ERR,ERROR,*)
724 
725  !Argument variables
726  TYPE(boundary_conditions_type), POINTER :: BOUNDARY_CONDITIONS
727  INTEGER(INTG), INTENT(OUT) :: ERR
728  TYPE(varying_string), INTENT(OUT) :: ERROR
729  !Local Variables
730 
731  enters("BOUNDARY_CONDITIONS_DESTROY",err,error,*999)
732 
733  IF(ASSOCIATED(boundary_conditions)) THEN
734  CALL boundary_conditions_finalise(boundary_conditions,err,error,*999)
735  ELSE
736  CALL flagerror("Boundary conditions is not associated.",err,error,*999)
737  ENDIF
738 
739  exits("BOUNDARY_CONDITIONS_DESTROY")
740  RETURN
741 999 errorsexits("BOUNDARY_CONDITIONS_DESTROY",err,error)
742  RETURN 1
743 
744  END SUBROUTINE boundary_conditions_destroy
745 
746  !
747  !================================================================================================================================
748  !
749 
751  SUBROUTINE boundary_conditions_finalise(BOUNDARY_CONDITIONS,ERR,ERROR,*)
752 
753  !Argument variables
754  TYPE(boundary_conditions_type), POINTER :: BOUNDARY_CONDITIONS
755  INTEGER(INTG), INTENT(OUT) :: ERR
756  TYPE(varying_string), INTENT(OUT) :: ERROR
757  !Local Variables
758  INTEGER(INTG) :: variable_idx
759 
760  enters("BOUNDARY_CONDITIONS_FINALISE",err,error,*999)
761 
762  IF(ASSOCIATED(boundary_conditions)) THEN
763  IF(ALLOCATED(boundary_conditions%BOUNDARY_CONDITIONS_VARIABLES)) THEN
764  DO variable_idx=1,boundary_conditions%NUMBER_OF_BOUNDARY_CONDITIONS_VARIABLES
765  IF(ASSOCIATED(boundary_conditions%BOUNDARY_CONDITIONS_VARIABLES(variable_idx)%PTR)) THEN
766  CALL boundary_conditions_variable_finalise(boundary_conditions%BOUNDARY_CONDITIONS_VARIABLES(variable_idx)%PTR, &
767  & err,error,*999)
768  ELSE
769  CALL flagerror("Boundary conditions variable number "//trim(number_to_vstring(variable_idx,"*",err,error))// &
770  & " is not associated",err,error,*999)
771  ENDIF
772  ENDDO !variable_idx
773  NULLIFY(boundary_conditions%SOLVER_EQUATIONS%SOLVER%SOLVER_EQUATIONS)
774  !BOUNDARY_CONDITIONS%SOLVER_EQUATIONS%SOLVER_EQUATIONS_FINISHED = .FALSE.
775  !BOUNDARY_CONDITIONS%SOLVER_EQUATIONS%SOLVER_MAPPING%SOLVER_MAPPING_FINISHED = .FALSE.
776  DEALLOCATE(boundary_conditions%BOUNDARY_CONDITIONS_VARIABLES)
777  ENDIF
778  DEALLOCATE(boundary_conditions)
779  ENDIF
780 
781  exits("BOUNDARY_CONDITIONS_FINALISE")
782  RETURN
783 999 errorsexits("BOUNDARY_CONDITIONS_FINALISE",err,error)
784  RETURN 1
785  END SUBROUTINE boundary_conditions_finalise
786 
787  !
788  !================================================================================================================================
789  !
790 
792  SUBROUTINE boundary_conditions_initialise(SOLVER_EQUATIONS,ERR,ERROR,*)
793 
794  !Argument variables
795  TYPE(solver_equations_type), POINTER :: SOLVER_EQUATIONS
796  INTEGER(INTG), INTENT(OUT) :: ERR
797  TYPE(varying_string), INTENT(OUT) :: ERROR
798  !Local Variables
799  INTEGER(INTG) :: DUMMY_ERR,variable_idx,variable_type,equations_set_idx,interface_condition_idx
800  TYPE(equations_type), POINTER :: EQUATIONS
801  TYPE(equations_set_type), POINTER :: EQUATIONS_SET
802  TYPE(equations_mapping_type), POINTER :: EQUATIONS_MAPPING
803  TYPE(equations_mapping_dynamic_type), POINTER :: DYNAMIC_MAPPING
804  TYPE(equations_mapping_linear_type), POINTER :: LINEAR_MAPPING
805  TYPE(equations_mapping_nonlinear_type), POINTER :: NONLINEAR_MAPPING
806  TYPE(equations_mapping_rhs_type), POINTER :: RHS_MAPPING
807  TYPE(interface_condition_type), POINTER :: INTERFACE_CONDITION
808  TYPE(interface_equations_type), POINTER :: INTERFACE_EQUATIONS
809  TYPE(interface_mapping_type), POINTER :: INTERFACE_MAPPING
810  TYPE(interface_mapping_rhs_type), POINTER :: INTERFACE_RHS_MAPPING
811  TYPE(varying_string) :: DUMMY_ERROR,LOCAL_ERROR
812 
813  enters("BOUNDARY_CONDITIONS_INITIALISE",err,error,*998)
814 
815  IF(ASSOCIATED(solver_equations)) THEN
816  IF(ASSOCIATED(solver_equations%BOUNDARY_CONDITIONS)) THEN
817  CALL flagerror("Boundary conditions is already associated for these solver equations.",err,error,*998)
818  ELSE
819  IF(ASSOCIATED(solver_equations%SOLVER_MAPPING)) THEN
820  ALLOCATE(solver_equations%BOUNDARY_CONDITIONS,stat=err)
821  IF(err/=0) CALL flagerror("Could not allocate boundary conditions.",err,error,*999)
822  solver_equations%BOUNDARY_CONDITIONS%BOUNDARY_CONDITIONS_FINISHED=.false.
823  solver_equations%BOUNDARY_CONDITIONS%NUMBER_OF_BOUNDARY_CONDITIONS_VARIABLES=0
824  solver_equations%BOUNDARY_CONDITIONS%SOLVER_EQUATIONS=>solver_equations
825  solver_equations%BOUNDARY_CONDITIONS%neumannMatrixSparsity=boundary_condition_sparse_matrices
826  DO equations_set_idx=1,solver_equations%SOLVER_MAPPING%NUMBER_OF_EQUATIONS_SETS
827  equations_set=>solver_equations%SOLVER_MAPPING%EQUATIONS_SETS(equations_set_idx)%PTR
828  IF(ASSOCIATED(equations_set)) THEN
829  equations=>equations_set%EQUATIONS
830  IF(ASSOCIATED(equations)) THEN
831  IF(equations%EQUATIONS_FINISHED) THEN
832  equations_mapping=>equations%EQUATIONS_MAPPING
833  IF(ASSOCIATED(equations_mapping)) THEN
834  IF(equations_mapping%EQUATIONS_MAPPING_FINISHED) THEN
835  equations_set%BOUNDARY_CONDITIONS=>solver_equations%BOUNDARY_CONDITIONS
836  SELECT CASE(equations%TIME_DEPENDENCE)
837  CASE(equations_static,equations_quasistatic)
838  SELECT CASE(equations%LINEARITY)
839  CASE(equations_linear,equations_nonlinear_bcs)
840  linear_mapping=>equations_mapping%LINEAR_MAPPING
841  IF(ASSOCIATED(linear_mapping)) THEN
842  DO variable_idx=1,linear_mapping%NUMBER_OF_LINEAR_MATRIX_VARIABLES
843  variable_type=linear_mapping%LINEAR_MATRIX_VARIABLE_TYPES(variable_idx)
844  IF(linear_mapping%VAR_TO_EQUATIONS_MATRICES_MAPS(variable_type)%NUMBER_OF_EQUATIONS_MATRICES>0) THEN
845  CALL boundary_conditions_variable_initialise(solver_equations%BOUNDARY_CONDITIONS, &
846  & linear_mapping%VAR_TO_EQUATIONS_MATRICES_MAPS(variable_type)%VARIABLE,err,error,*999)
847  ENDIF
848  ENDDO !variable_idx
849  ELSE
850  CALL flagerror("Equations mapping linear mapping is not associated.",err,error,*999)
851  ENDIF
852  rhs_mapping=>equations_mapping%RHS_MAPPING
853  IF(ASSOCIATED(rhs_mapping)) THEN
854  CALL boundary_conditions_variable_initialise(solver_equations%BOUNDARY_CONDITIONS, &
855  & rhs_mapping%RHS_VARIABLE,err,error,*999)
856  ENDIF
857  CASE(equations_nonlinear)
858  nonlinear_mapping=>equations_mapping%NONLINEAR_MAPPING
859  IF(ASSOCIATED(nonlinear_mapping)) THEN
860  DO variable_idx=1,nonlinear_mapping%NUMBER_OF_RESIDUAL_VARIABLES
861  CALL boundary_conditions_variable_initialise(solver_equations%BOUNDARY_CONDITIONS, &
862  & nonlinear_mapping%RESIDUAL_VARIABLES(variable_idx)%PTR,err,error,*999)
863  ENDDO
864  ELSE
865  CALL flagerror("Equations mapping nonlinear mapping is not associated.",err,error,*999)
866  ENDIF
867  rhs_mapping=>equations_mapping%RHS_MAPPING
868  IF(ASSOCIATED(rhs_mapping)) THEN
869  CALL boundary_conditions_variable_initialise(solver_equations%BOUNDARY_CONDITIONS, &
870  & rhs_mapping%RHS_VARIABLE,err,error,*999)
871  ELSE
872  CALL flagerror("Equations mapping RHS mapping is not associated.",err,error,*999)
873  ENDIF
874  CASE DEFAULT
875  local_error="The equations linearity type of "//trim(number_to_vstring(equations%LINEARITY,"*", &
876  & err,error))//" is invalid."
877  CALL flagerror(local_error,err,error,*999)
878  END SELECT
879  CASE(equations_first_order_dynamic,equations_second_order_dynamic)
880  SELECT CASE(equations%LINEARITY)
881  CASE(equations_linear,equations_nonlinear_bcs)
882  dynamic_mapping=>equations_mapping%DYNAMIC_MAPPING
883  IF(ASSOCIATED(dynamic_mapping)) THEN
884  CALL boundary_conditions_variable_initialise(solver_equations%BOUNDARY_CONDITIONS, &
885  & dynamic_mapping%DYNAMIC_VARIABLE,err,error,*999)
886  ELSE
887  CALL flagerror("Equations mapping dynamic mapping is not associated.",err,error,*999)
888  ENDIF
889  rhs_mapping=>equations_mapping%RHS_MAPPING
890  IF(ASSOCIATED(rhs_mapping)) THEN
891  CALL boundary_conditions_variable_initialise(solver_equations%BOUNDARY_CONDITIONS, &
892  & rhs_mapping%RHS_VARIABLE,err,error,*999)
893  ELSE
894  CALL flagerror("Equations mapping RHS mapping is not associated.",err,error,*999)
895  ENDIF
896  CASE(equations_nonlinear)
897  dynamic_mapping=>equations_mapping%DYNAMIC_MAPPING
898  IF(ASSOCIATED(dynamic_mapping)) THEN
899  CALL boundary_conditions_variable_initialise(solver_equations%BOUNDARY_CONDITIONS, &
900  & dynamic_mapping%DYNAMIC_VARIABLE,err,error,*999)
901  ELSE
902  CALL flagerror("Equations mapping dynamic mapping is not associated.",err,error,*999)
903  ENDIF
904  rhs_mapping=>equations_mapping%RHS_MAPPING
905  IF(ASSOCIATED(rhs_mapping)) THEN
906  CALL boundary_conditions_variable_initialise(solver_equations%BOUNDARY_CONDITIONS, &
907  & rhs_mapping%RHS_VARIABLE,err,error,*999)
908  ELSE
909  CALL flagerror("Equations mapping RHS mapping is not associated.",err,error,*999)
910  ENDIF
911  CASE DEFAULT
912  local_error="The equations linearity type of "//trim(number_to_vstring(equations%LINEARITY,"*", &
913  & err,error))//" is invalid."
914  CALL flagerror(local_error,err,error,*999)
915  END SELECT
916  CASE DEFAULT
917  local_error="The equations time dependence type of "// &
918  & trim(number_to_vstring(equations%TIME_DEPENDENCE,"*",err,error))//" is invalid."
919  CALL flagerror(local_error,err,error,*999)
920  END SELECT
921  ELSE
922  CALL flagerror("Equations mapping has not been finished.",err,error,*998)
923  ENDIF
924  ELSE
925  CALL flagerror("Equations equations mapping is not associated.",err,error,*998)
926  ENDIF
927  ELSE
928  CALL flagerror("Equations has not been finished.",err,error,*998)
929  ENDIF
930  ELSE
931  CALL flagerror("Equations set equations is not associated.",err,error,*998)
932  ENDIF
933  ELSE
934  CALL flagerror("Equations set is not associated.",err,error,*998)
935  ENDIF
936  ENDDO !equations_set_idx
937  DO interface_condition_idx=1,solver_equations%SOLVER_MAPPING%NUMBER_OF_INTERFACE_CONDITIONS
938  interface_condition=>solver_equations%SOLVER_MAPPING%INTERFACE_CONDITIONS(interface_condition_idx)%PTR
939  IF(ASSOCIATED(interface_condition)) THEN
940  interface_equations=>interface_condition%INTERFACE_EQUATIONS
941  IF(ASSOCIATED(interface_equations)) THEN
942  IF(interface_equations%INTERFACE_EQUATIONS_FINISHED) THEN
943  interface_mapping=>interface_equations%INTERFACE_MAPPING
944  IF(ASSOCIATED(interface_mapping)) THEN
945  IF(interface_mapping%INTERFACE_MAPPING_FINISHED) THEN
946  interface_condition%BOUNDARY_CONDITIONS=>solver_equations%BOUNDARY_CONDITIONS
947  !Only linear interface equations implemented at the moment
948  SELECT CASE(interface_equations%TIME_DEPENDENCE)
949  CASE(interface_condition_static,interface_condition_quasistatic)
950  SELECT CASE(interface_equations%LINEARITY)
951  CASE(interface_condition_linear)
952  interface_mapping=>interface_equations%INTERFACE_MAPPING
953  IF(ASSOCIATED(interface_mapping)) THEN
954  variable_type=interface_mapping%LAGRANGE_VARIABLE_TYPE
955  IF(interface_mapping%NUMBER_OF_INTERFACE_MATRICES>0) THEN
956  CALL boundary_conditions_variable_initialise(solver_equations%BOUNDARY_CONDITIONS, &
957  & interface_mapping%LAGRANGE_VARIABLE,err,error,*999)
958  ENDIF
959  ELSE
960  CALL flagerror("Interface mapping mapping is not associated.",err,error,*999)
961  ENDIF
962  interface_rhs_mapping=>interface_mapping%RHS_MAPPING
963  IF(ASSOCIATED(interface_rhs_mapping)) THEN
964  CALL boundary_conditions_variable_initialise(solver_equations%BOUNDARY_CONDITIONS, &
965  & interface_rhs_mapping%RHS_VARIABLE,err,error,*999)
966  ENDIF
967  CASE DEFAULT
968  local_error="The equations linearity type of "//trim(number_to_vstring(equations%LINEARITY,"*", &
969  & err,error))//" is invalid."
970  CALL flagerror(local_error,err,error,*999)
971  END SELECT
972  CASE DEFAULT
973  local_error="The equations time dependence type of "// &
974  & trim(number_to_vstring(equations%TIME_DEPENDENCE,"*",err,error))//" is invalid."
975  CALL flagerror(local_error,err,error,*999)
976  END SELECT
977  ELSE
978  CALL flagerror("Interface mapping has not been finished.",err,error,*998)
979  ENDIF
980  ELSE
981  CALL flagerror("Interface mapping is not associated.",err,error,*998)
982  ENDIF
983  ELSE
984  CALL flagerror("Interface equations has not been finished.",err,error,*998)
985  ENDIF
986  ELSE
987  CALL flagerror("Interface equations is not associated.",err,error,*998)
988  ENDIF
989  ELSE
990  CALL flagerror("Interface condition not associated.",err,error,*998)
991  ENDIF
992  ENDDO !interface_condition_idx
993  ELSE
994  CALL flagerror("Solver equations solver mapping is not associated.",err,error,*998)
995  ENDIF
996  ENDIF
997  ELSE
998  CALL flagerror("Solver equations is not associated",err,error,*998)
999  ENDIF
1000 
1001  exits("BOUNDARY_CONDITIONS_INITIALISE")
1002  RETURN
1003 999 CALL boundary_conditions_finalise(solver_equations%BOUNDARY_CONDITIONS,dummy_err,dummy_error,*998)
1004 998 errorsexits("BOUNDARY_CONDITIONS_INITIALISE",err,error)
1005  RETURN 1
1006 
1007  END SUBROUTINE boundary_conditions_initialise
1008 
1009  !
1010  !================================================================================================================================
1011  !
1012 
1014  SUBROUTINE boundary_conditions_add_constant(BOUNDARY_CONDITIONS,FIELD,VARIABLE_TYPE,COMPONENT_NUMBER,CONDITION,VALUE,ERR,ERROR,*)
1015 
1016  !Argument variables
1017  TYPE(boundary_conditions_type), POINTER :: BOUNDARY_CONDITIONS
1018  TYPE(field_type), POINTER :: FIELD
1019  INTEGER(INTG), INTENT(IN) :: VARIABLE_TYPE
1020  INTEGER(INTG), INTENT(IN) :: COMPONENT_NUMBER
1021  INTEGER(INTG), INTENT(IN) :: CONDITION
1022  REAL(DP), INTENT(IN) :: VALUE
1023  INTEGER(INTG), INTENT(OUT) :: ERR
1024  TYPE(varying_string), INTENT(OUT) :: ERROR
1025  !Local Variables
1026  INTEGER(INTG) :: local_ny,global_ny
1027  TYPE(boundary_conditions_variable_type), POINTER :: BOUNDARY_CONDITIONS_VARIABLE
1028  TYPE(field_variable_type), POINTER :: DEPENDENT_FIELD_VARIABLE
1029  TYPE(varying_string) :: LOCAL_ERROR
1030 
1031  enters("BOUNDARY_CONDITIONS_ADD_CONSTANT",err,error,*999)
1032 
1033  NULLIFY(boundary_conditions_variable)
1034  NULLIFY(dependent_field_variable)
1035 
1036  !Note: This routine is for constant interpolation
1037  IF(ASSOCIATED(boundary_conditions)) THEN
1038  IF(boundary_conditions%BOUNDARY_CONDITIONS_FINISHED) THEN
1039  CALL flagerror("Boundary conditions have been finished.",err,error,*999)
1040  ELSE
1041  IF(ASSOCIATED(field)) THEN
1042  CALL field_component_dof_get_constant(field,variable_type,component_number,local_ny,global_ny, &
1043  & err,error,*999)
1044  CALL field_variable_get(field,variable_type,dependent_field_variable,err,error,*999)
1045  CALL boundary_conditions_variable_get(boundary_conditions,dependent_field_variable,boundary_conditions_variable, &
1046  & err,error,*999)
1047  IF(ASSOCIATED(boundary_conditions_variable)) THEN
1048  CALL boundaryconditions_checkinterpolationtype(condition,field,variable_type,component_number,err,error,*999)
1049  CALL boundary_conditions_add_local_dof(boundary_conditions,field,variable_type, &
1050  & local_ny,condition,VALUE,err,error,*999)
1051  ELSE
1052  local_error="The boundary conditions for variable type "//trim(number_to_vstring(variable_type,"*",err,error))// &
1053  & " has not been created."
1054  CALL flagerror(local_error,err,error,*999)
1055  ENDIF
1056  ELSE
1057  CALL flagerror("The dependent field is not associated.",err,error,*999)
1058  ENDIF
1059  ENDIF
1060  ELSE
1061  CALL flagerror("Boundary conditions is not associated.",err,error,*999)
1062  ENDIF
1063 
1064  exits("BOUNDARY_CONDITION_ADD_CONSTANT")
1065  RETURN
1066 999 errorsexits("BOUNDARY_CONDITION_ADD_CONSTANT",err,error)
1067  RETURN 1
1068  END SUBROUTINE boundary_conditions_add_constant
1069 
1070  !
1071  !================================================================================================================================
1072  !
1073 
1075  SUBROUTINE boundary_conditions_set_constant(BOUNDARY_CONDITIONS,FIELD,VARIABLE_TYPE,COMPONENT_NUMBER,CONDITION,VALUE,ERR,ERROR,*)
1076 
1077  !Argument variables
1078  TYPE(boundary_conditions_type), POINTER :: BOUNDARY_CONDITIONS
1079  TYPE(field_type), POINTER :: FIELD
1080  INTEGER(INTG), INTENT(IN) :: VARIABLE_TYPE
1081  INTEGER(INTG), INTENT(IN) :: COMPONENT_NUMBER
1082  INTEGER(INTG), INTENT(IN) :: CONDITION
1083  REAL(DP), INTENT(IN) :: VALUE
1084  INTEGER(INTG), INTENT(OUT) :: ERR
1085  TYPE(varying_string), INTENT(OUT) :: ERROR
1086  !Local Variables
1087  INTEGER(INTG) :: local_ny,global_ny
1088  TYPE(boundary_conditions_variable_type), POINTER :: BOUNDARY_CONDITIONS_VARIABLE
1089  TYPE(field_variable_type), POINTER :: FIELD_VARIABLE
1090  TYPE(varying_string) :: LOCAL_ERROR
1091 
1092  enters("BOUNDARY_CONDITIONS_SET_CONSTANT",err,error,*999)
1093 
1094  !Note: This routine is for constant interpolation
1095  IF(ASSOCIATED(boundary_conditions)) THEN
1096  IF(boundary_conditions%BOUNDARY_CONDITIONS_FINISHED) THEN
1097  CALL flagerror("Boundary conditions have been finished.",err,error,*999)
1098  ELSE
1099  IF(ASSOCIATED(field)) THEN
1100  CALL field_component_dof_get_constant(field,variable_type,component_number,local_ny,global_ny, &
1101  & err,error,*999)
1102  CALL field_variable_get(field,variable_type,field_variable,err,error,*999)
1103  CALL boundary_conditions_variable_get(boundary_conditions,field_variable,boundary_conditions_variable,err,error,*999)
1104  IF(ASSOCIATED(boundary_conditions_variable)) THEN
1105  CALL boundaryconditions_checkinterpolationtype(condition,field,variable_type,component_number,err,error,*999)
1106  CALL boundary_conditions_set_local_dof(boundary_conditions,field,variable_type, &
1107  & local_ny,condition,VALUE,err,error,*999)
1108  ELSE
1109  local_error="The boundary conditions for variable type "//trim(number_to_vstring(variable_type,"*",err,error))// &
1110  & " has not been created."
1111  CALL flagerror(local_error,err,error,*999)
1112  ENDIF
1113  ELSE
1114  CALL flagerror("The dependent field is not associated.",err,error,*999)
1115  ENDIF
1116  ENDIF
1117  ELSE
1118  CALL flagerror("Boundary conditions is not associated.",err,error,*999)
1119  ENDIF
1120 
1121  exits("BOUNDARY_CONDITION_SET_CONSTANT")
1122  RETURN
1123 999 errorsexits("BOUNDARY_CONDITION_SET_CONSTANT",err,error)
1124  RETURN 1
1125  END SUBROUTINE boundary_conditions_set_constant
1126 
1127  !
1128  !================================================================================================================================
1129  !
1130 
1132  SUBROUTINE boundary_conditions_add_local_dof1(BOUNDARY_CONDITIONS,FIELD,VARIABLE_TYPE,DOF_INDEX,CONDITION,VALUE,ERR,ERROR,*)
1133 
1134  !Argument variables
1135  TYPE(boundary_conditions_type), POINTER :: BOUNDARY_CONDITIONS
1136  TYPE(field_type), POINTER :: FIELD
1137  INTEGER(INTG), INTENT(IN) :: VARIABLE_TYPE
1138  INTEGER(INTG), INTENT(IN) :: DOF_INDEX
1139  INTEGER(INTG), INTENT(IN) :: CONDITION
1140  REAL(DP), INTENT(IN) :: VALUE
1141  INTEGER(INTG), INTENT(OUT) :: ERR
1142  TYPE(varying_string), INTENT(OUT) :: ERROR
1143  !Local Variables
1144 
1145  enters("BOUNDARY_CONDITIONS_ADD_LOCAL_DOF1",err,error,*999)
1146 
1147  CALL boundary_conditions_add_local_dofs(boundary_conditions,field,variable_type,(/dof_index/),(/condition/),(/VALUE/), &
1148  & err,error,*999)
1149 
1150  exits("BOUNDARY_CONDITIONS_ADD_LOCAL_DOF1")
1151  RETURN
1152 999 errorsexits("BOUNDARY_CONDITIONS_ADD_LOCAL_DOF1",err,error)
1153  RETURN 1
1154  END SUBROUTINE boundary_conditions_add_local_dof1
1155 
1156  !
1157  !================================================================================================================================
1158  !
1159 
1161  SUBROUTINE boundary_conditions_add_local_dofs(BOUNDARY_CONDITIONS,FIELD,VARIABLE_TYPE,DOF_INDICES,CONDITIONS,VALUES,ERR,ERROR,*)
1162 
1163  !Argument variables
1164  TYPE(boundary_conditions_type), POINTER :: BOUNDARY_CONDITIONS
1165  TYPE(field_type), POINTER :: FIELD
1166  INTEGER(INTG), INTENT(IN) :: VARIABLE_TYPE
1167  INTEGER(INTG), INTENT(IN) :: DOF_INDICES(:)
1168  INTEGER(INTG), INTENT(IN) :: CONDITIONS(:)
1169  REAL(DP), INTENT(IN) :: VALUES(:)
1170  INTEGER(INTG), INTENT(OUT) :: ERR
1171  TYPE(varying_string), INTENT(OUT) :: ERROR
1172  !Local Variables
1173  INTEGER(INTG) :: i,global_ny,local_ny
1174  REAL(DP) :: INITIAL_VALUE
1175  TYPE(boundary_conditions_variable_type), POINTER :: BOUNDARY_CONDITIONS_VARIABLE
1176  TYPE(domain_mapping_type), POINTER :: DOMAIN_MAPPING
1177  TYPE(field_variable_type), POINTER :: DEPENDENT_VARIABLE
1178  TYPE(varying_string) :: LOCAL_ERROR
1179 
1180  enters("BOUNDARY_CONDITIONS_ADD_LOCAL_DOFS",err,error,*999)
1181  NULLIFY(dependent_variable)
1182 
1183  IF(ASSOCIATED(boundary_conditions)) THEN
1184  IF(boundary_conditions%BOUNDARY_CONDITIONS_FINISHED) THEN
1185  CALL flagerror("Boundary conditions have been finished.",err,error,*999)
1186  ELSE
1187  IF(ASSOCIATED(field)) THEN
1188  NULLIFY(dependent_variable)
1189  CALL field_variable_get(field,variable_type,dependent_variable,err,error,*999)
1190  IF(ASSOCIATED(dependent_variable)) THEN
1191  domain_mapping=>dependent_variable%DOMAIN_MAPPING
1192  IF(ASSOCIATED(domain_mapping)) THEN
1193  CALL boundary_conditions_variable_get(boundary_conditions,dependent_variable,boundary_conditions_variable, &
1194  & err,error,*999)
1195  IF(ASSOCIATED(boundary_conditions_variable)) THEN
1196  IF(SIZE(dof_indices,1)==SIZE(conditions,1)) THEN
1197  IF(SIZE(dof_indices,1)==SIZE(values,1)) THEN
1198  DO i=1,SIZE(dof_indices,1)
1199  local_ny=dof_indices(i)
1200  IF(local_ny>=1.AND.local_ny<=domain_mapping%NUMBER_OF_LOCAL) THEN
1201  global_ny=domain_mapping%LOCAL_TO_GLOBAL_MAP(local_ny)
1202  ! Set boundary condition and dof type, and make sure parameter sets are created
1203  CALL boundaryconditions_setconditiontype(boundary_conditions_variable,global_ny,conditions(i), &
1204  & err,error,*999)
1205  ! Update field sets by adding boundary condition values
1206  SELECT CASE(conditions(i))
1207  CASE(boundary_condition_free)
1208  ! No field update
1209  CASE(boundary_condition_fixed)
1210  CALL field_parameter_set_add_local_dof(field,variable_type,field_values_set_type, &
1211  & local_ny,values(i),err,error,*999)
1212  CASE(boundary_condition_fixed_inlet)
1213  CALL field_parameter_set_add_local_dof(field,variable_type,field_values_set_type, &
1214  & local_ny,values(i),err,error,*999)
1215  CASE(boundary_condition_fixed_outlet)
1216  CALL field_parameter_set_add_local_dof(field,variable_type,field_values_set_type, &
1217  & local_ny,values(i),err,error,*999)
1218  CASE(boundary_condition_fixed_wall)
1219  CALL field_parameter_set_add_local_dof(field,variable_type,field_values_set_type, &
1220  & local_ny,values(i),err,error,*999)
1221  CASE(boundary_condition_moved_wall)
1222  CALL field_parameter_set_add_local_dof(field,variable_type,field_values_set_type, &
1223  & local_ny,values(i),err,error,*999)
1224  CASE(boundary_condition_free_wall)
1225  CALL field_parameter_set_add_local_dof(field,variable_type,field_values_set_type, &
1226  & local_ny,values(i),err,error,*999)
1227  CASE(boundary_condition_moved_wall_incremented)
1228  CALL field_parameter_set_add_local_dof(field,variable_type,field_values_set_type,local_ny,values(i), &
1229  & err,error,*999)
1230  CALL field_parameter_set_add_local_dof(field,variable_type,field_boundary_conditions_set_type, &
1231  & local_ny,values(i),err,error,*999)
1232  CASE(boundary_condition_fixed_incremented)
1233  ! For increment loops, we need to set the full BC parameter set value by
1234  ! getting the current value from the values parameter set
1235  CALL field_parameter_set_get_local_dof(field,variable_type,field_values_set_type, &
1236  & local_ny,initial_value,err,error,*999)
1237  CALL field_parameter_set_update_local_dof(field,variable_type,field_boundary_conditions_set_type, &
1238  & local_ny,initial_value+values(i),err,error,*999)
1239  CASE(boundary_condition_pressure)
1240  CALL field_parameter_set_add_local_dof(field,variable_type,field_pressure_values_set_type, &
1241  & local_ny,values(i),err,error,*999)
1242  CASE(boundary_condition_pressure_incremented)
1243  ! For pressure incremented, adding to the values_set parameter value doesn't make sense,
1244  ! so just increment the value in the pressure values parameter set
1245  CALL field_parameter_set_add_local_dof(field,variable_type,field_pressure_values_set_type, &
1246  & local_ny,values(i),err,error,*999)
1247  CASE(boundary_condition_correction_mass_increase)
1248  ! No field update
1249  CASE(boundary_condition_impermeable_wall)
1250  CALL field_parameter_set_add_local_dof(field,variable_type,field_impermeable_flag_values_set_type, &
1251  & local_ny,values(i),err,error,*999)
1252  CASE(boundary_condition_neumann_point)
1253  ! Point value is stored in boundary conditions field set, and is then integrated to
1254  ! get the RHS variable value
1255  CALL field_parameter_set_add_local_dof(field,variable_type,field_boundary_conditions_set_type, &
1256  & local_ny,values(i),err,error,*999)
1257  CASE(boundary_condition_neumann_point_incremented)
1258  CALL field_parameter_set_add_local_dof(field,variable_type,field_boundary_conditions_set_type, &
1259  & local_ny,values(i),err,error,*999)
1260  CASE(boundary_condition_neumann_integrated,boundary_condition_neumann_integrated_only)
1261  ! For integrated Neumann condition, integration is already done, so set the RHS
1262  ! dof value directly
1263  CALL field_parameter_set_add_local_dof(field,variable_type,field_values_set_type, &
1264  & local_ny,values(i),err,error,*999)
1265  CASE(boundary_condition_fixed_fitted,boundary_condition_fixed_nonreflecting, &
1266  & boundary_condition_fixed_cellml,boundary_condition_fixed_stree)
1267  CALL field_parameter_set_add_local_dof(field,variable_type,field_values_set_type, &
1268  & local_ny,values(i),err,error,*999)
1269  CASE DEFAULT
1270  local_error="The specified boundary condition type for dof index "// &
1271  & trim(number_to_vstring(i,"*",err,error))//" of "// &
1272  & trim(number_to_vstring(conditions(i),"*",err,error))//" is invalid."
1273  CALL flagerror(local_error,err,error,*999)
1274  END SELECT
1275  ELSE
1276  local_error="The local dof of "//&
1277  & trim(number_to_vstring(local_ny,"*",err,error))//" at dof index "// &
1278  & trim(number_to_vstring(i,"*",err,error))// &
1279  & " is invalid. The dof should be between 1 and "// &
1280  & trim(number_to_vstring(domain_mapping%NUMBER_OF_LOCAL,"*",err,error))//"."
1281  CALL flagerror(local_error,err,error,*999)
1282  ENDIF
1283  ENDDO !i
1284  ELSE
1285  local_error="The size of the dof indices array ("// &
1286  & trim(number_to_vstring(SIZE(dof_indices,1),"*",err,error))// &
1287  & ") does not match the size of the values array ("// &
1288  & trim(number_to_vstring(SIZE(values,1),"*",err,error))//")."
1289  CALL flagerror(local_error,err,error,*999)
1290  ENDIF
1291  ELSE
1292  local_error="The size of the dof indices array ("// &
1293  & trim(number_to_vstring(SIZE(dof_indices,1),"*",err,error))// &
1294  & ") does not match the size of the fixed conditions array ("// &
1295  & trim(number_to_vstring(SIZE(conditions,1),"*",err,error))//")."
1296  CALL flagerror(local_error,err,error,*999)
1297  ENDIF
1298  ELSE
1299  CALL flagerror("Boundary conditions variable is not associated.",err,error,*999)
1300  ENDIF
1301  ELSE
1302  CALL flagerror("The dependent field variable domain mapping is not associated.",err,error,*999)
1303  ENDIF
1304  ELSE
1305  CALL flagerror("The dependent field variable is not associated.",err,error,*999)
1306  ENDIF
1307  ELSE
1308  CALL flagerror("The dependent field is not associated..",err,error,*999)
1309  ENDIF
1310  ENDIF
1311  ELSE
1312  CALL flagerror("Boundary conditions is not associated.",err,error,*999)
1313  ENDIF
1314 
1315  exits("BOUNDARY_CONDITIONS_ADD_LOCAL_DOFS")
1316  RETURN
1317 999 errorsexits("BOUNDARY_CONDITIONS_ADD_LOCAL_DOFS",err,error)
1318  RETURN 1
1319  END SUBROUTINE boundary_conditions_add_local_dofs
1320 
1321  !
1322  !================================================================================================================================
1323  !
1324 
1326  SUBROUTINE boundary_conditions_set_local_dof1(BOUNDARY_CONDITIONS,FIELD,VARIABLE_TYPE,DOF_INDEX,CONDITION,VALUE,ERR,ERROR,*)
1327 
1328  !Argument variables
1329  TYPE(boundary_conditions_type), POINTER :: BOUNDARY_CONDITIONS
1330  TYPE(field_type), POINTER :: FIELD
1331  INTEGER(INTG), INTENT(IN) :: VARIABLE_TYPE
1332  INTEGER(INTG), INTENT(IN) :: DOF_INDEX
1333  INTEGER(INTG), INTENT(IN) :: CONDITION
1334  REAL(DP), INTENT(IN) :: VALUE
1335  INTEGER(INTG), INTENT(OUT) :: ERR
1336  TYPE(varying_string), INTENT(OUT) :: ERROR
1337  !Local Variables
1338 
1339  enters("BOUNDARY_CONDITIONS_SET_LOCAL_DOF1",err,error,*999)
1340  CALL boundary_conditions_set_local_dofs(boundary_conditions,field,variable_type,(/dof_index/),(/condition/),(/VALUE/), &
1341  & err,error,*999)
1342 
1343  exits("BOUNDARY_CONDITIONS_SET_LOCAL_DOF1")
1344  RETURN
1345 999 errorsexits("BOUNDARY_CONDITIONS_SET_LOCAL_DOF1",err,error)
1346  RETURN 1
1347  END SUBROUTINE boundary_conditions_set_local_dof1
1348 
1349  !
1350  !================================================================================================================================
1351  !
1352 
1354  SUBROUTINE boundary_conditions_set_local_dofs(BOUNDARY_CONDITIONS,FIELD,VARIABLE_TYPE,DOF_INDICES,CONDITIONS,VALUES,ERR,ERROR,*)
1355 
1356  !Argument variables
1357  TYPE(boundary_conditions_type), POINTER :: BOUNDARY_CONDITIONS
1358  TYPE(field_type), POINTER :: FIELD
1359  INTEGER(INTG), INTENT(IN) :: VARIABLE_TYPE
1360  INTEGER(INTG), INTENT(IN) :: DOF_INDICES(:)
1361  INTEGER(INTG), INTENT(IN) :: CONDITIONS(:)
1362  REAL(DP), INTENT(IN) :: VALUES(:)
1363  INTEGER(INTG), INTENT(OUT) :: ERR
1364  TYPE(varying_string), INTENT(OUT) :: ERROR
1365  !Local Variables
1366  INTEGER(INTG) :: i,global_ny,local_ny
1367  TYPE(boundary_conditions_variable_type), POINTER :: BOUNDARY_CONDITIONS_VARIABLE
1368  TYPE(domain_mapping_type), POINTER :: DOMAIN_MAPPING
1369  TYPE(field_variable_type), POINTER :: DEPENDENT_VARIABLE
1370  TYPE(varying_string) :: LOCAL_ERROR
1371 
1372  enters("BOUNDARY_CONDITIONS_SET_LOCAL_DOFS",err,error,*999)
1373 
1374  IF(ASSOCIATED(boundary_conditions)) THEN
1375  IF(boundary_conditions%BOUNDARY_CONDITIONS_FINISHED) THEN
1376  CALL flagerror("Boundary conditions have been finished.",err,error,*999)
1377  ELSE
1378  IF(ASSOCIATED(field)) THEN
1379  NULLIFY(dependent_variable)
1380  CALL field_variable_get(field,variable_type,dependent_variable,err,error,*999)
1381  IF(ASSOCIATED(dependent_variable)) THEN
1382  domain_mapping=>dependent_variable%DOMAIN_MAPPING
1383  IF(ASSOCIATED(domain_mapping)) THEN
1384  CALL boundary_conditions_variable_get(boundary_conditions,dependent_variable,boundary_conditions_variable, &
1385  & err,error,*999)
1386  IF(ASSOCIATED(boundary_conditions_variable)) THEN
1387  IF(SIZE(dof_indices,1)==SIZE(conditions,1)) THEN
1388  IF(SIZE(dof_indices,1)==SIZE(values,1)) THEN
1389  DO i=1,SIZE(dof_indices,1)
1390  local_ny=dof_indices(i)
1391  IF(local_ny>=1.AND.local_ny<=domain_mapping%NUMBER_OF_LOCAL) THEN
1392  global_ny=domain_mapping%LOCAL_TO_GLOBAL_MAP(local_ny)
1393  ! Set boundary condition and dof type
1394  CALL boundaryconditions_setconditiontype(boundary_conditions_variable,global_ny,conditions(i), &
1395  & err,error,*999)
1396  ! Update field sets with boundary condition value
1397 
1398  SELECT CASE(conditions(i))
1399  CASE(boundary_condition_free)
1400  ! No field update
1401  CASE(boundary_condition_fixed)
1402  CALL field_parameter_set_update_local_dof(field,variable_type,field_values_set_type, &
1403  & local_ny,values(i),err,error,*999)
1404  CASE(boundary_condition_fixed_inlet)
1405  CALL field_parameter_set_update_local_dof(field,variable_type,field_values_set_type, &
1406  & local_ny,values(i),err,error,*999)
1407  CASE(boundary_condition_fixed_outlet)
1408  CALL field_parameter_set_update_local_dof(field,variable_type,field_values_set_type, &
1409  & local_ny,values(i),err,error,*999)
1410  CASE(boundary_condition_fixed_wall)
1411  CALL field_parameter_set_update_local_dof(field,variable_type,field_values_set_type, &
1412  & local_ny,values(i),err,error,*999)
1413  CASE(boundary_condition_moved_wall)
1414  CALL field_parameter_set_update_local_dof(field,variable_type,field_values_set_type, &
1415  & local_ny,values(i),err,error,*999)
1416  CASE(boundary_condition_free_wall)
1417  CALL field_parameter_set_update_local_dof(field,variable_type,field_values_set_type, &
1418  & local_ny,values(i),err,error,*999)
1419  CASE(boundary_condition_moved_wall_incremented)
1420  CALL field_parameter_set_update_local_dof(field,variable_type,field_values_set_type,local_ny,values(i), &
1421  & err,error,*999)
1422  CALL field_parameter_set_update_local_dof(field,variable_type,field_boundary_conditions_set_type, &
1423  & local_ny,values(i),err,error,*999)
1424  CASE(boundary_condition_fixed_incremented) !For load increment loops
1425  CALL field_parameter_set_update_local_dof(field,variable_type,field_boundary_conditions_set_type, &
1426  & local_ny,values(i),err,error,*999)
1427  CASE(boundary_condition_pressure)
1428  CALL field_parameter_set_update_local_dof(field,variable_type,field_pressure_values_set_type, &
1429  & local_ny,values(i),err,error,*999)
1430  CASE(boundary_condition_pressure_incremented)
1431  CALL field_parameter_set_update_local_dof(field,variable_type,field_pressure_values_set_type, &
1432  & local_ny,values(i),err,error,*999)
1433  CASE(boundary_condition_correction_mass_increase)
1434  ! No field update
1435  CASE(boundary_condition_neumann_point)
1436  CALL field_parameter_set_update_local_dof(field,variable_type,field_boundary_conditions_set_type, &
1437  & local_ny,values(i),err,error,*999)
1438  CASE(boundary_condition_neumann_point_incremented)
1439  CALL field_parameter_set_update_local_dof(field,variable_type,field_boundary_conditions_set_type, &
1440  & local_ny,values(i),err,error,*999)
1441  CASE(boundary_condition_neumann_integrated,boundary_condition_neumann_integrated_only)
1442  CALL field_parameter_set_update_local_dof(field,variable_type,field_values_set_type, &
1443  & local_ny,values(i),err,error,*999)
1444  CASE(boundary_condition_impermeable_wall)
1445  CALL field_parameter_set_update_local_dof(field,variable_type,field_impermeable_flag_values_set_type, &
1446  & local_ny,values(i),err,error,*999)
1447  CASE(boundary_condition_fixed_fitted,boundary_condition_fixed_nonreflecting, &
1448  & boundary_condition_fixed_cellml,boundary_condition_fixed_stree)
1449  CALL field_parameter_set_update_local_dof(field,variable_type,field_values_set_type, &
1450  & local_ny,values(i),err,error,*999)
1451  CASE DEFAULT
1452  local_error="The specified boundary condition type for dof index "// &
1453  & trim(number_to_vstring(i,"*",err,error))//" of "// &
1454  & trim(number_to_vstring(conditions(i),"*",err,error))//" is invalid."
1455  CALL flagerror(local_error,err,error,*999)
1456  END SELECT
1457  ELSE
1458  local_error="The local dof of "//&
1459  & trim(number_to_vstring(local_ny,"*",err,error))//" at dof index "// &
1460  & trim(number_to_vstring(i,"*",err,error))// &
1461  & " is invalid. The dof should be between 1 and "// &
1462  & trim(number_to_vstring(domain_mapping%NUMBER_OF_LOCAL,"*",err,error))//"."
1463  CALL flagerror(local_error,err,error,*999)
1464  ENDIF
1465  ENDDO !i
1466  ELSE
1467  local_error="The size of the dof indices array ("// &
1468  & trim(number_to_vstring(SIZE(dof_indices,1),"*",err,error))// &
1469  & ") does not match the size of the values array ("// &
1470  & trim(number_to_vstring(SIZE(values,1),"*",err,error))//")."
1471  CALL flagerror(local_error,err,error,*999)
1472  ENDIF
1473  ELSE
1474  local_error="The size of the dof indices array ("// &
1475  & trim(number_to_vstring(SIZE(dof_indices,1),"*",err,error))// &
1476  & ") does not match the size of the fixed conditions array ("// &
1477  & trim(number_to_vstring(SIZE(conditions,1),"*",err,error))//")."
1478  CALL flagerror(local_error,err,error,*999)
1479  ENDIF
1480  ELSE
1481  CALL flagerror("Boundary conditions variable is not associated.",err,error,*999)
1482  ENDIF
1483  ELSE
1484  CALL flagerror("The dependent field variable domain mapping is not associated.",err,error,*999)
1485  ENDIF
1486  ELSE
1487  CALL flagerror("The dependent field variable is not associated.",err,error,*999)
1488  ENDIF
1489  ELSE
1490  CALL flagerror("The dependent field is not associated.",err,error,*999)
1491  ENDIF
1492  ENDIF
1493  ELSE
1494  CALL flagerror("Boundary conditions is not associated.",err,error,*999)
1495  ENDIF
1496 
1497  exits("BOUNDARY_CONDITIONS_SET_LOCAL_DOFS")
1498  RETURN
1499 999 errorsexits("BOUNDARY_CONDITIONS_SET_LOCAL_DOFS",err,error)
1500  RETURN 1
1501  END SUBROUTINE boundary_conditions_set_local_dofs
1502 
1503  !
1504  !================================================================================================================================
1505  !
1506 
1509  SUBROUTINE boundaryconditions_setconditiontype(boundaryConditionsVariable,globalDof,condition,err,error,*)
1510 
1511  !Argument variables
1512  TYPE(boundary_conditions_variable_type), POINTER :: boundaryConditionsVariable
1513  INTEGER(INTG), INTENT(IN) :: globalDof
1514  INTEGER(INTG), INTENT(IN) :: condition
1515  INTEGER(INTG), INTENT(OUT) :: err
1516  TYPE(varying_string), INTENT(OUT) :: error
1517  !Local variables
1518  INTEGER(INTG) :: dofType, previousCondition, previousDof
1519 
1520  enters("BoundaryConditions_SetConditionType",err,error,*999)
1521 
1522  ! We won't do much checking here as this is only used internally and everything has been checked for
1523  ! association already
1524  ! Don't need to make sure field values set type is available as this will always be there, but need
1525  ! to make sure any other parameter sets required are.
1526  SELECT CASE(condition)
1527  CASE(boundary_condition_free)
1528  doftype=boundary_condition_dof_free
1529  CASE(boundary_condition_fixed)
1530  doftype=boundary_condition_dof_fixed
1531  CASE(boundary_condition_linear_constraint)
1532  doftype=boundary_condition_dof_constrained
1533  CASE(boundary_condition_fixed_inlet)
1534  doftype=boundary_condition_dof_fixed
1535  CASE(boundary_condition_fixed_outlet)
1536  doftype=boundary_condition_dof_fixed
1537  CASE(boundary_condition_fixed_wall)
1538  doftype=boundary_condition_dof_fixed
1539  CASE(boundary_condition_moved_wall)
1540  doftype=boundary_condition_dof_fixed
1541  CASE(boundary_condition_free_wall)
1542  doftype=boundary_condition_dof_free
1543  CASE(boundary_condition_moved_wall_incremented)
1544  doftype=boundary_condition_dof_fixed
1545  CALL field_parametersetensurecreated(boundaryconditionsvariable%VARIABLE%FIELD,boundaryconditionsvariable%VARIABLE_TYPE, &
1546  & field_boundary_conditions_set_type,err,error,*999)
1547  boundaryconditionsvariable%parameterSetRequired(field_boundary_conditions_set_type)=.true.
1548  CASE(boundary_condition_fixed_incremented) !For load increment loops
1549  doftype=boundary_condition_dof_fixed
1550  CALL field_parametersetensurecreated(boundaryconditionsvariable%VARIABLE%FIELD,boundaryconditionsvariable%VARIABLE_TYPE, &
1551  & field_boundary_conditions_set_type,err,error,*999)
1552  boundaryconditionsvariable%parameterSetRequired(field_boundary_conditions_set_type)=.true.
1553  CASE(boundary_condition_pressure)
1554  ! Pressure boundary conditions leave the RHS dof as free, as the Neumann terms
1555  ! are calculated in finite elasticity routines when calculating the element residual
1556  doftype=boundary_condition_dof_free
1557  CALL field_parametersetensurecreated(boundaryconditionsvariable%VARIABLE%FIELD,boundaryconditionsvariable%VARIABLE_TYPE, &
1558  & field_pressure_values_set_type,err,error,*999)
1559  boundaryconditionsvariable%parameterSetRequired(field_pressure_values_set_type)=.true.
1560  CASE(boundary_condition_pressure_incremented)
1561  doftype=boundary_condition_dof_free
1562  CALL field_parametersetensurecreated(boundaryconditionsvariable%VARIABLE%FIELD,boundaryconditionsvariable%VARIABLE_TYPE, &
1563  & field_pressure_values_set_type,err,error,*999)
1564  boundaryconditionsvariable%parameterSetRequired(field_pressure_values_set_type)=.true.
1565  CALL field_parametersetensurecreated(boundaryconditionsvariable%VARIABLE%FIELD,boundaryconditionsvariable%VARIABLE_TYPE, &
1566  & field_previous_pressure_set_type,err,error,*999)
1567  boundaryconditionsvariable%parameterSetRequired(field_previous_pressure_set_type)=.true.
1568  CASE(boundary_condition_correction_mass_increase)
1569  doftype=boundary_condition_dof_fixed
1570  CASE(boundary_condition_impermeable_wall)
1571  doftype=boundary_condition_dof_free
1572  CALL field_parametersetensurecreated(boundaryconditionsvariable%VARIABLE%FIELD,boundaryconditionsvariable%VARIABLE_TYPE, &
1573  & field_impermeable_flag_values_set_type,err,error,*999)
1574  boundaryconditionsvariable%parameterSetRequired(field_impermeable_flag_values_set_type)=.true.
1575  CASE(boundary_condition_neumann_point,boundary_condition_neumann_point_incremented)
1576  doftype=boundary_condition_dof_fixed
1577  CALL field_parametersetensurecreated(boundaryconditionsvariable%VARIABLE%FIELD,boundaryconditionsvariable%VARIABLE_TYPE, &
1578  & field_boundary_conditions_set_type,err,error,*999)
1579  CALL field_parametersetensurecreated(boundaryconditionsvariable%VARIABLE%FIELD,boundaryconditionsvariable%VARIABLE_TYPE, &
1580  & field_integrated_neumann_set_type,err,error,*999)
1581  boundaryconditionsvariable%parameterSetRequired(field_boundary_conditions_set_type)=.true.
1582  boundaryconditionsvariable%parameterSetRequired(field_integrated_neumann_set_type)=.true.
1583  CASE(boundary_condition_neumann_integrated,boundary_condition_neumann_integrated_only)
1584  doftype=boundary_condition_dof_fixed
1585  CASE(boundary_condition_fixed_fitted,boundary_condition_fixed_nonreflecting,boundary_condition_fixed_cellml, &
1586  & boundary_condition_fixed_stree)
1587  doftype=boundary_condition_dof_fixed
1588  CASE DEFAULT
1589  CALL flagerror("The specified boundary condition type for dof number "// &
1590  & trim(number_to_vstring(globaldof,"*",err,error))//" of "// &
1591  & trim(number_to_vstring(condition,"*",err,error))//" is invalid.", &
1592  & err,error,*999)
1593  END SELECT
1594 
1595  !We have a valid boundary condition type
1596  !Update condition type counts
1597  previouscondition=boundaryconditionsvariable%CONDITION_TYPES(globaldof)
1598  IF(previouscondition/=condition) THEN
1599  ! DOF_COUNTS array doesn't include a count for BOUNDARY_CONDITION_FREE, which equals zero
1600  IF(previouscondition/=boundary_condition_free) THEN
1601  boundaryconditionsvariable%DOF_COUNTS(previouscondition)= &
1602  & boundaryconditionsvariable%DOF_COUNTS(previouscondition)-1
1603  END IF
1604  IF(condition/=boundary_condition_free) THEN
1605  boundaryconditionsvariable%DOF_COUNTS(condition)= &
1606  & boundaryconditionsvariable%DOF_COUNTS(condition)+1
1607  END IF
1608  END IF
1609  !Update Dirichlet DOF count
1610  previousdof=boundaryconditionsvariable%DOF_TYPES(globaldof)
1611  IF(doftype==boundary_condition_dof_fixed.AND.previousdof/=boundary_condition_dof_fixed) THEN
1612  boundaryconditionsvariable%NUMBER_OF_DIRICHLET_CONDITIONS= &
1613  & boundaryconditionsvariable%NUMBER_OF_DIRICHLET_CONDITIONS+1
1614  ELSE IF(doftype/=boundary_condition_dof_fixed.AND.previousdof==boundary_condition_dof_fixed) THEN
1615  boundaryconditionsvariable%NUMBER_OF_DIRICHLET_CONDITIONS= &
1616  & boundaryconditionsvariable%NUMBER_OF_DIRICHLET_CONDITIONS-1
1617  END IF
1618 
1619  !Set the boundary condition type and DOF type
1620  boundaryconditionsvariable%CONDITION_TYPES(globaldof)=condition
1621  boundaryconditionsvariable%DOF_TYPES(globaldof)=doftype
1622  IF(diagnostics1) THEN
1623  CALL write_string(diagnostic_output_type,"Boundary Condition Being Set",err,error,*999)
1624  CALL write_string_value(diagnostic_output_type,"global dof = ", globaldof,err,error,*999)
1625  CALL write_string_value(diagnostic_output_type,"Variable Type = ", &
1626  & boundaryconditionsvariable%VARIABLE_TYPE,err,error,*999)
1627  CALL write_string_value(diagnostic_output_type,"New Condition = ", &
1628  & condition,err,error,*999)
1629  CALL write_string_value(diagnostic_output_type,"dof type = ", &
1630  & doftype,err,error,*999)
1631  ENDIF
1632  exits("BoundaryConditions_SetConditionType")
1633  RETURN
1634 999 errorsexits("BoundaryConditions_SetConditionType",err,error)
1635  RETURN 1
1636  END SUBROUTINE boundaryconditions_setconditiontype
1637 
1638  !
1639  !================================================================================================================================
1640  !
1641 
1643  SUBROUTINE boundary_conditions_add_element(BOUNDARY_CONDITIONS,FIELD,VARIABLE_TYPE,USER_ELEMENT_NUMBER,COMPONENT_NUMBER, &
1644  & condition,VALUE,err,error,*)
1645 
1646  !Argument variables
1647  TYPE(boundary_conditions_type), POINTER :: BOUNDARY_CONDITIONS
1648  TYPE(field_type), POINTER :: FIELD
1649  INTEGER(INTG), INTENT(IN) :: VARIABLE_TYPE
1650  INTEGER(INTG), INTENT(IN) :: USER_ELEMENT_NUMBER
1651  INTEGER(INTG), INTENT(IN) :: COMPONENT_NUMBER
1652  INTEGER(INTG), INTENT(IN) :: CONDITION
1653  REAL(DP), INTENT(IN) :: VALUE
1654  INTEGER(INTG), INTENT(OUT) :: ERR
1655  TYPE(varying_string), INTENT(OUT) :: ERROR
1656  !Local Variables
1657  INTEGER(INTG) :: local_ny,global_ny
1658  TYPE(boundary_conditions_variable_type), POINTER :: BOUNDARY_CONDITIONS_VARIABLE
1659  TYPE(field_variable_type), POINTER :: FIELD_VARIABLE
1660  TYPE(varying_string) :: LOCAL_ERROR
1661 
1662  enters("BOUNDARY_CONDITIONS_ADD_ELEMENT",err,error,*999)
1663 
1664  !Note: this routine is for element based interpolation
1665  IF(ASSOCIATED(boundary_conditions)) THEN
1666  IF(boundary_conditions%BOUNDARY_CONDITIONS_FINISHED) THEN
1667  CALL flagerror("Boundary conditions have been finished.",err,error,*999)
1668  ELSE
1669  IF(ASSOCIATED(field)) THEN
1670  CALL field_component_dof_get_user_element(field,variable_type,user_element_number,component_number, &
1671  & local_ny,global_ny,err,error,*999)
1672  NULLIFY(field_variable)
1673  NULLIFY(boundary_conditions_variable)
1674  CALL field_variable_get(field,variable_type,field_variable,err,error,*999)
1675  IF(ASSOCIATED(field_variable)) THEN
1676  CALL boundary_conditions_variable_get(boundary_conditions,field_variable,boundary_conditions_variable,err,error,*999)
1677  IF(ASSOCIATED(boundary_conditions_variable)) THEN
1678  CALL boundaryconditions_checkinterpolationtype(condition,field,variable_type,component_number,err,error,*999)
1679  CALL boundary_conditions_add_local_dof(boundary_conditions,field,variable_type, &
1680  & local_ny,condition,VALUE,err,error,*999)
1681  ELSE
1682  local_error="The boundary conditions for variable type "//trim(number_to_vstring(variable_type,"*",err,error))// &
1683  & " has not been created."
1684  CALL flagerror(local_error,err,error,*999)
1685  ENDIF
1686  ELSE
1687  CALL flagerror("The dependent field variable is not associated.",err,error,*999)
1688  ENDIF
1689  ELSE
1690  CALL flagerror("The dependent field is not associated.",err,error,*999)
1691  ENDIF
1692  ENDIF
1693  ELSE
1694  CALL flagerror("Boundary conditions is not associated.",err,error,*999)
1695  ENDIF
1696 
1697  exits("BOUNDARY_CONDITION_ADD_ELEMENT")
1698  RETURN
1699 999 errorsexits("BOUNDARY_CONDITION_ADD_ELEMENT",err,error)
1700  RETURN 1
1701  END SUBROUTINE boundary_conditions_add_element
1702 
1703  !
1704  !================================================================================================================================
1705  !
1706 
1708  SUBROUTINE boundaryconditions_checkinterpolationtype(condition,field,variableType,componentNumber,err,error,*)
1709 
1710  ! Argument variables
1711  INTEGER(INTG), INTENT(IN) :: condition
1712  TYPE(field_type), POINTER :: field
1713  INTEGER(INTG), INTENT(IN) :: variableType
1714  INTEGER(INTG), INTENT(IN) :: componentNumber
1715  INTEGER(INTG), INTENT(OUT) :: err
1716  TYPE(varying_string), INTENT(OUT) :: error
1717  ! Local variables
1718  INTEGER(INTG) :: interpolationType
1719  LOGICAL :: validCondition
1720 
1721  enters("BoundaryConditions_CheckInterpolationType",err,error,*999)
1722 
1723  CALL field_component_interpolation_get(field,variabletype,componentnumber,interpolationtype,err,error,*999)
1724 
1725  validcondition=.true.
1726  SELECT CASE(condition)
1727  CASE(boundary_condition_free, &
1728  & boundary_condition_fixed, &
1729  & boundary_condition_fixed_incremented)
1730  ! Valid for all interpolation types
1731  CASE(boundary_condition_fixed_inlet, &
1732  & boundary_condition_fixed_outlet)
1733  IF(interpolationtype/=field_node_based_interpolation) THEN
1734  validcondition=.false.
1735  END IF
1736  CASE(boundary_condition_fixed_wall, &
1737  & boundary_condition_moved_wall, &
1738  & boundary_condition_free_wall, &
1739  & boundary_condition_moved_wall_incremented)
1740  IF(interpolationtype/=field_node_based_interpolation) THEN
1741  validcondition=.false.
1742  END IF
1743  CASE(boundary_condition_pressure, &
1744  & boundary_condition_pressure_incremented)
1745  IF(interpolationtype/=field_node_based_interpolation) THEN
1746  validcondition=.false.
1747  END IF
1748  CASE(boundary_condition_correction_mass_increase)
1749  IF(interpolationtype/=field_node_based_interpolation) THEN
1750  validcondition=.false.
1751  END IF
1752  CASE(boundary_condition_impermeable_wall)
1753  IF(interpolationtype/=field_node_based_interpolation) THEN
1754  validcondition=.false.
1755  END IF
1756  CASE(boundary_condition_neumann_point, &
1757  & boundary_condition_neumann_integrated, &
1758  & boundary_condition_neumann_integrated_only, &
1759  & boundary_condition_neumann_point_incremented)
1760  IF(interpolationtype/=field_node_based_interpolation) THEN
1761  validcondition=.false.
1762  END IF
1763  CASE(boundary_condition_fixed_fitted,boundary_condition_fixed_nonreflecting,boundary_condition_fixed_cellml, &
1764  & boundary_condition_fixed_stree)
1765  IF(interpolationtype/=field_node_based_interpolation) THEN
1766  validcondition=.false.
1767  END IF
1768  CASE DEFAULT
1769  CALL flagerror("The specified boundary condition type of "// &
1770  & trim(number_to_vstring(condition,"*",err,error))//" is invalid.", &
1771  & err,error,*999)
1772  END SELECT
1773  IF(.NOT.validcondition) THEN
1774  CALL flagerror("The specified boundary condition type of "// &
1775  & trim(number_to_vstring(condition,"*",err,error))//" is not valid for the field component "// &
1776  & "interpolation type of "//trim(number_to_vstring(interpolationtype,"*",err,error))//".", &
1777  & err,error,*999)
1778  END IF
1779 
1780  exits("BoundaryConditions_CheckInterpolationType")
1781  RETURN
1782 999 errorsexits("BoundaryConditions_CheckInterpolationType",err,error)
1783  RETURN 1
1784  END SUBROUTINE boundaryconditions_checkinterpolationtype
1785 
1786  !
1787  !================================================================================================================================
1788  !
1789 
1791  SUBROUTINE boundaryconditions_checkequations(boundaryConditionsVariable,err,error,*)
1792 
1793  ! Argument variables
1794  TYPE(boundary_conditions_variable_type), POINTER :: boundaryConditionsVariable
1795  INTEGER(INTG), INTENT(OUT) :: err
1796  type(varying_string), intent(out) :: error
1797  ! Local variables
1798  INTEGER(INTG) :: boundaryConditionType,equationsSetIdx,specificationSize
1799  TYPE(solver_equations_type), POINTER :: solverEquations
1800  TYPE(solver_mapping_type), POINTER :: solverMapping
1801  TYPE(equations_set_type), POINTER :: equationsSet
1802  LOGICAL :: validEquationsSetFound
1803 
1804  enters("BoundaryConditions_CheckEquations",err,error,*999)
1805 
1806  !Get and check pointers we need
1807  solverequations=>boundaryconditionsvariable%BOUNDARY_CONDITIONS%SOLVER_EQUATIONS
1808  IF(.NOT.ASSOCIATED(solverequations)) THEN
1809  CALL flagerror("Boundary conditions solver equations are not associated.",err,error,*999)
1810  END IF
1811  solvermapping=>solverequations%SOLVER_MAPPING
1812  IF(.NOT.ASSOCIATED(solvermapping)) THEN
1813  CALL flagerror("Solver equations solver mapping is not associated.",err,error,*999)
1814  END IF
1815 
1816  DO boundaryconditiontype=1,max_boundary_condition_number
1817  !Check if any DOFs have been set to this BC type
1818  IF(boundaryconditionsvariable%DOF_COUNTS(boundaryconditiontype)>0) THEN
1819  validequationssetfound=.false.
1820  DO equationssetidx=1,solvermapping%NUMBER_OF_EQUATIONS_SETS
1821  equationsset=>solvermapping%EQUATIONS_SETS(equationssetidx)%PTR
1822  IF(.NOT.ASSOCIATED(equationsset)) THEN
1823  CALL flagerror("Solver equations equations set is not associated.",err,error,*999)
1824  END IF
1825  IF(.NOT.ALLOCATED(equationsset%specification)) THEN
1826  CALL flagerror("Equations set specification is not allocated.",err,error,*999)
1827  END IF
1828  specificationsize=SIZE(equationsset%specification,1)
1829 
1830  SELECT CASE(boundaryconditiontype)
1831  CASE(boundary_condition_free)
1832  ! Valid for any equations set
1833  validequationssetfound=.true.
1834  CASE(boundary_condition_fixed)
1835  validequationssetfound=.true.
1836  CASE(boundary_condition_linear_constraint)
1837  validequationssetfound=.true.
1838  CASE(boundary_condition_fixed_inlet, &
1839  & boundary_condition_fixed_outlet)
1840  IF(specificationsize>=2) THEN
1841  IF(equationsset%specification(1)==equations_set_fluid_mechanics_class.AND. &
1842  & (equationsset%specification(2)==equations_set_stokes_equation_type.OR. &
1843  & equationsset%specification(2)==equations_set_characteristic_equation_type.OR. &
1844  & equationsset%specification(2)==equations_set_navier_stokes_equation_type)) THEN
1845  validequationssetfound=.true.
1846  END IF
1847  END IF
1848  CASE(boundary_condition_fixed_wall,boundary_condition_moved_wall, &
1849  & boundary_condition_moved_wall_incremented,boundary_condition_free_wall)
1850  IF(specificationsize>=2) THEN
1851  IF(equationsset%specification(1)==equations_set_fluid_mechanics_class.AND. &
1852  & (equationsset%specification(2)==equations_set_stokes_equation_type.OR. &
1853  & equationsset%specification(2)==equations_set_navier_stokes_equation_type.OR. &
1854  & equationsset%specification(2)==equations_set_darcy_equation_type)) THEN
1855  validequationssetfound=.true.
1856  ELSE IF(specificationsize==3) THEN
1857  IF(equationsset%specification(1)==equations_set_classical_field_class.AND. &
1858  & equationsset%specification(2)==equations_set_laplace_equation_type.AND. &
1859  & equationsset%specification(3)==equations_set_moving_mesh_laplace_subtype) THEN
1860  validequationssetfound=.true.
1861  END IF
1862  END IF
1863  END IF
1864  CASE(boundary_condition_fixed_incremented)
1865  validequationssetfound=.true.
1866  CASE(boundary_condition_pressure, &
1867  & boundary_condition_pressure_incremented)
1868  IF(specificationsize>=2) THEN
1869  IF(equationsset%specification(1)==equations_set_elasticity_class.AND. &
1870  & equationsset%specification(2)==equations_set_finite_elasticity_type) THEN
1871  validequationssetfound=.true.
1872  ELSE IF(equationsset%specification(1)==equations_set_fluid_mechanics_class .AND. &
1873  & equationsset%specification(2)==equations_set_navier_stokes_equation_type) THEN
1874  validequationssetfound=.true.
1875  END IF
1876  ENDIF
1877  CASE(boundary_condition_correction_mass_increase)
1878  !Not actually used anywhere? So keep it as invalid, although maybe it should be removed?
1879  validequationssetfound=.false.
1880  CASE(boundary_condition_impermeable_wall)
1881  IF(specificationsize>=3) THEN
1882  IF(equationsset%specification(1)==equations_set_elasticity_class.AND. &
1883  & equationsset%specification(2)==equations_set_finite_elasticity_type.AND. &
1884  & (equationsset%specification(3)==equations_set_incompressible_finite_elasticity_darcy_subtype.OR. &
1885  & equationsset%specification(3)==equations_set_elasticity_darcy_inria_model_subtype.OR. &
1886  & equationsset%specification(3)==equations_set_incompressible_elasticity_driven_darcy_subtype)) THEN
1887  validequationssetfound=.true.
1888  END IF
1889  END IF
1890  CASE(boundary_condition_neumann_point,boundary_condition_neumann_point_incremented)
1891  validequationssetfound=.true.
1892  CASE(boundary_condition_neumann_integrated,boundary_condition_neumann_integrated_only)
1893  validequationssetfound=.true.
1894  CASE(boundary_condition_fixed_fitted)
1895  IF(equationsset%specification(1)==equations_set_fluid_mechanics_class.AND. &
1896  & (equationsset%specification(2)==equations_set_stokes_equation_type.OR. &
1897  & equationsset%specification(2)==equations_set_navier_stokes_equation_type)) THEN
1898  validequationssetfound=.true.
1899  END IF
1900  CASE(boundary_condition_fixed_nonreflecting,boundary_condition_fixed_cellml,boundary_condition_fixed_stree)
1901  IF(equationsset%specification(1)==equations_set_fluid_mechanics_class.AND. &
1902  & (equationsset%specification(2)==equations_set_characteristic_equation_type.OR. &
1903  & equationsset%specification(2)==equations_set_navier_stokes_equation_type)) THEN
1904  validequationssetfound=.true.
1905  END IF
1906  CASE DEFAULT
1907  CALL flagerror("The specified boundary condition type of "// &
1908  & trim(number_to_vstring(boundaryconditiontype,"*",err,error))// &
1909  & " is invalid.",err,error,*999)
1910  END SELECT
1911  END DO
1912 
1913  IF(.NOT.validequationssetfound) THEN
1914  CALL flagerror("The specified boundary condition type of "// &
1915  & trim(number_to_vstring(boundaryconditiontype,"*",err,error))// &
1916  & " is invalid for the equations sets in the solver equations.",err,error,*999)
1917  END IF
1918  END IF
1919  END DO
1920 
1921  exits("BoundaryConditions_CheckEquations")
1922  RETURN
1923 999 errorsexits("BoundaryConditions_CheckEquations",err,error)
1924  RETURN 1
1925  END SUBROUTINE boundaryconditions_checkequations
1926 
1927  !
1928  !================================================================================================================================
1929  !
1930 
1932  SUBROUTINE boundary_conditions_set_element(BOUNDARY_CONDITIONS,FIELD,VARIABLE_TYPE,USER_ELEMENT_NUMBER,COMPONENT_NUMBER, &
1933  & condition,VALUE,err,error,*)
1934 
1935  !Argument variables
1936  TYPE(boundary_conditions_type), POINTER :: BOUNDARY_CONDITIONS
1937  TYPE(field_type), POINTER :: FIELD
1938  INTEGER(INTG), INTENT(IN) :: VARIABLE_TYPE
1939  INTEGER(INTG), INTENT(IN) :: USER_ELEMENT_NUMBER
1940  INTEGER(INTG), INTENT(IN) :: COMPONENT_NUMBER
1941  INTEGER(INTG), INTENT(IN) :: CONDITION
1942  REAL(DP), INTENT(IN) :: VALUE
1943  INTEGER(INTG), INTENT(OUT) :: ERR
1944  TYPE(varying_string), INTENT(OUT) :: ERROR
1945  !Local Variables
1946  INTEGER(INTG) :: local_ny,global_ny
1947  TYPE(boundary_conditions_variable_type), POINTER :: BOUNDARY_CONDITIONS_VARIABLE
1948  TYPE(field_variable_type), POINTER :: FIELD_VARIABLE
1949  TYPE(varying_string) :: LOCAL_ERROR
1950 
1951  enters("BOUNDARY_CONDITIONS_SET_ELEMENT",err,error,*999)
1952 
1953  !Note: this routine is for element based interpolation
1954  IF(ASSOCIATED(boundary_conditions)) THEN
1955  IF(boundary_conditions%BOUNDARY_CONDITIONS_FINISHED) THEN
1956  CALL flagerror("Boundary conditions have been finished.",err,error,*999)
1957  ELSE
1958  IF(ASSOCIATED(field)) THEN
1959  CALL field_component_dof_get_user_element(field,variable_type,user_element_number,component_number, &
1960  & local_ny,global_ny,err,error,*999)
1961  NULLIFY(field_variable)
1962  NULLIFY(boundary_conditions_variable)
1963  CALL field_variable_get(field,variable_type,field_variable,err,error,*999)
1964  IF(ASSOCIATED(field_variable)) THEN
1965  CALL boundary_conditions_variable_get(boundary_conditions,field_variable,boundary_conditions_variable,err,error,*999)
1966  IF(ASSOCIATED(boundary_conditions_variable)) THEN
1967  CALL boundaryconditions_checkinterpolationtype(condition,field,variable_type,component_number,err,error,*999)
1968  CALL boundary_conditions_set_local_dof(boundary_conditions,field,variable_type, &
1969  & local_ny,condition,VALUE,err,error,*999)
1970  ELSE
1971  local_error="The boundary conditions for variable type "//trim(number_to_vstring(variable_type,"*",err,error))// &
1972  & " has not been created."
1973  CALL flagerror(local_error,err,error,*999)
1974  ENDIF
1975  ELSE
1976  CALL flagerror("The dependent field variable is not associated.",err,error,*999)
1977  ENDIF
1978  ELSE
1979  CALL flagerror("The dependent field is not associated.",err,error,*999)
1980  ENDIF
1981  ENDIF
1982  ELSE
1983  CALL flagerror("Boundary conditions is not associated.",err,error,*999)
1984  ENDIF
1985 
1986  exits("BOUNDARY_CONDITION_SET_ELEMENT")
1987  RETURN
1988 999 errorsexits("BOUNDARY_CONDITION_SET_ELEMENT",err,error)
1989  RETURN 1
1990  END SUBROUTINE boundary_conditions_set_element
1991 
1992  !
1993  !================================================================================================================================
1994  !
1995 
1997  SUBROUTINE boundary_conditions_add_node(BOUNDARY_CONDITIONS,FIELD,VARIABLE_TYPE,VERSION_NUMBER,DERIVATIVE_NUMBER, &
1998  & user_node_number,component_number,condition,VALUE,err,error,*)
1999 
2000  !Argument variables
2001  TYPE(boundary_conditions_type), POINTER :: BOUNDARY_CONDITIONS
2002  TYPE(field_type), POINTER :: FIELD
2003  INTEGER(INTG), INTENT(IN) :: VARIABLE_TYPE
2004  INTEGER(INTG), INTENT(IN) :: VERSION_NUMBER
2005  INTEGER(INTG), INTENT(IN) :: DERIVATIVE_NUMBER
2006  INTEGER(INTG), INTENT(IN) :: USER_NODE_NUMBER
2007  INTEGER(INTG), INTENT(IN) :: COMPONENT_NUMBER
2008  INTEGER(INTG), INTENT(IN) :: CONDITION
2009  REAL(DP), INTENT(IN) :: VALUE
2010  INTEGER(INTG), INTENT(OUT) :: ERR
2011  TYPE(varying_string), INTENT(OUT) :: ERROR
2012  !Local Variables
2013  INTEGER(INTG) :: local_ny,global_ny
2014  TYPE(boundary_conditions_variable_type), POINTER :: BOUNDARY_CONDITIONS_VARIABLE
2015  TYPE(field_variable_type), POINTER :: FIELD_VARIABLE
2016  TYPE(varying_string) :: LOCAL_ERROR
2017 
2018  enters("BOUNDARY_CONDITIONS_ADD_NODE",err,error,*999)
2019 
2020  NULLIFY(field_variable)
2021  NULLIFY(boundary_conditions_variable)
2022 
2023  IF(ASSOCIATED(boundary_conditions)) THEN
2024  IF(boundary_conditions%BOUNDARY_CONDITIONS_FINISHED) THEN
2025  CALL flagerror("Boundary conditions have been finished.",err,error,*999)
2026  ELSE
2027  IF(ASSOCIATED(field)) THEN
2028  CALL field_component_dof_get_user_node(field,variable_type,version_number,derivative_number, &
2029  & user_node_number,component_number,local_ny,global_ny,err,error,*999)
2030  CALL field_variable_get(field,variable_type,field_variable,err,error,*999)
2031  IF(ASSOCIATED(field_variable)) THEN
2032  CALL boundary_conditions_variable_get(boundary_conditions,field_variable,boundary_conditions_variable,err,error,*999)
2033  IF(ASSOCIATED(boundary_conditions_variable)) THEN
2034  CALL boundaryconditions_checkinterpolationtype(condition,field,variable_type,component_number,err,error,*999)
2035  CALL boundary_conditions_add_local_dof(boundary_conditions,field,variable_type, &
2036  & local_ny,condition,VALUE,err,error,*999)
2037  ELSE
2038  local_error="The boundary conditions for variable type "//trim(number_to_vstring(variable_type,"*",err,error))// &
2039  & " has not been created."
2040  CALL flagerror(local_error,err,error,*999)
2041  ENDIF
2042  ELSE
2043  CALL flagerror("The dependent field variable is not associated",err,error,*999)
2044  ENDIF
2045  ELSE
2046  CALL flagerror("The dependent field is not associated.",err,error,*999)
2047  ENDIF
2048  ENDIF
2049  ELSE
2050  CALL flagerror("Boundary conditions is not associated.",err,error,*999)
2051  ENDIF
2052 
2053  exits("BOUNDARY_CONDITIONS_ADD_NODE")
2054  RETURN
2055 999 errorsexits("BOUNDARY_CONDITIONS_ADD_NODE",err,error)
2056  RETURN 1
2057  END SUBROUTINE boundary_conditions_add_node
2058 
2059  !
2060  !================================================================================================================================
2061  !
2062 
2064  SUBROUTINE boundaryconditions_neumanninitialise(boundaryConditionsVariable,err,error,*)
2065 
2066  !Argument variables
2067  TYPE(boundary_conditions_variable_type), POINTER :: boundaryConditionsVariable
2068  INTEGER(INTG), INTENT(OUT) :: err
2069  TYPE(varying_string), INTENT(OUT) :: error
2070  !Local Variables
2071  TYPE(boundaryconditionsneumanntype), POINTER :: boundaryConditionsNeumann
2072  INTEGER(INTG) :: numberOfValues,numberOfLocalDofs
2073  INTEGER(INTG) :: dummyErr
2074  TYPE(varying_string) :: dummyError
2075 
2076  enters("BoundaryConditions_NeumannInitialise",err,error,*998)
2077 
2078  IF(ASSOCIATED(boundaryconditionsvariable)) THEN
2079  numberofvalues=boundaryconditionsvariable%DOF_COUNTS(boundary_condition_neumann_point)+ &
2080  & boundaryconditionsvariable%DOF_COUNTS(boundary_condition_neumann_point_incremented)
2081  ALLOCATE(boundaryconditionsvariable%neumannBoundaryConditions,stat=err)
2082  IF(err/=0) CALL flagerror("Could not allocate Neumann Boundary Conditions",err,error,*998)
2083  boundaryconditionsneumann=>boundaryconditionsvariable%neumannBoundaryConditions
2084  IF(ASSOCIATED(boundaryconditionsneumann)) THEN
2085  NULLIFY(boundaryconditionsneumann%integrationMatrix)
2086  NULLIFY(boundaryconditionsneumann%pointValues)
2087  NULLIFY(boundaryconditionsneumann%pointDofMapping)
2088 
2089  numberoflocaldofs=boundaryconditionsvariable%VARIABLE%NUMBER_OF_DOFS
2090  ALLOCATE(boundaryconditionsneumann%setDofs(numberofvalues),stat=err)
2091  IF(err/=0) CALL flagerror("Could not allocate Neumann set DOFs.",err,error,*999)
2092  boundaryconditionsneumann%setDofs=0
2093  ELSE
2094  CALL flagerror("The boundary condition Neumann is not associated",err,error,*998)
2095  END IF
2096  ELSE
2097  CALL flagerror("Boundary conditions variable is not associated.",err,error,*998)
2098  END IF
2099 
2100  exits("BoundaryConditions_NeumannInitialise")
2101  RETURN
2102 999 CALL boundaryconditions_neumannfinalise(boundaryconditionsvariable,dummyerr,dummyerror,*998)
2103 998 errorsexits("BoundaryConditions_NeumannInitialise",err,error)
2104  RETURN 1
2105  END SUBROUTINE boundaryconditions_neumanninitialise
2106 
2107  !
2108  !================================================================================================================================
2109  !
2110 
2114  SUBROUTINE boundaryconditions_neumannmatricesinitialise(boundaryConditionsVariable,err,error,*)
2115 
2116  !Argument variables
2117  TYPE(boundary_conditions_variable_type), POINTER :: boundaryConditionsVariable
2118  INTEGER(INTG), INTENT(OUT) :: err
2119  TYPE(varying_string), INTENT(OUT) :: error
2120  !Local Variables
2121  TYPE(boundaryconditionsneumanntype), POINTER :: boundaryConditionsNeumann
2122  TYPE(field_variable_type), POINTER :: rhsVariable
2123  TYPE(domain_mapping_type), POINTER :: rowMapping, pointDofMapping
2124  TYPE(domain_topology_type), POINTER :: topology
2125  TYPE(domain_line_type), POINTER :: line
2126  TYPE(domain_face_type), POINTER :: face
2127  TYPE(list_type), POINTER :: columnIndicesList, rowColumnIndicesList
2128  INTEGER(INTG) :: myComputationalNodeNumber
2129  INTEGER(INTG) :: numberOfPointDofs, numberNonZeros, numberRowEntries, neumannConditionNumber, localNeumannConditionIdx
2130  INTEGER(INTG) :: neumannIdx, globalDof, localDof, localDofNyy, domainIdx, numberOfDomains, domainNumber, componentNumber
2131  INTEGER(INTG) :: nodeIdx, derivIdx, nodeNumber, versionNumber, derivativeNumber, columnNodeNumber, lineIdx, faceIdx, columnDof
2132  INTEGER(INTG), ALLOCATABLE :: rowIndices(:), columnIndices(:), localDofNumbers(:)
2133  REAL(DP) :: pointValue
2134  INTEGER(INTG) :: dummyErr
2135  TYPE(varying_string) :: dummyError
2136 
2137  enters("BoundaryConditions_NeumannMatricesInitialise",err,error,*999)
2138 
2139  IF(ASSOCIATED(boundaryconditionsvariable)) THEN
2140  rhsvariable=>boundaryconditionsvariable%variable
2141  IF(.NOT.ASSOCIATED(rhsvariable)) &
2142  & CALL flagerror("RHS boundary conditions variable field variable is not associated.",err,error,*999)
2143  numberofpointdofs=boundaryconditionsvariable%DOF_COUNTS(boundary_condition_neumann_point) + &
2144  & boundaryconditionsvariable%DOF_COUNTS(boundary_condition_neumann_point_incremented)
2145  boundaryconditionsneumann=>boundaryconditionsvariable%neumannBoundaryConditions
2146  IF(ASSOCIATED(boundaryconditionsneumann)) THEN
2147  ! For rows we can re-use the RHS variable row mapping
2148  rowmapping=>rhsvariable%DOMAIN_MAPPING
2149  IF(.NOT.ASSOCIATED(rowmapping)) &
2150  & CALL flagerror("RHS field variable mapping is not associated.",err,error,*998)
2151 
2152  ! Create a domain mapping for the Neumann point DOFs, required for the distributed matrix columns
2153  ALLOCATE(pointdofmapping,stat=err)
2154  IF(err/=0) CALL flagerror("Could not allocate Neumann DOF domain mapping.",err,error,*999)
2155  CALL domain_mappings_mapping_initialise(pointdofmapping,rowmapping%NUMBER_OF_DOMAINS,err,error,*999)
2156  boundaryconditionsneumann%pointDofMapping=>pointdofmapping
2157  ! Calculate global to local mapping for Neumann DOFs
2158  pointdofmapping%NUMBER_OF_GLOBAL=numberofpointdofs
2159  ALLOCATE(pointdofmapping%GLOBAL_TO_LOCAL_MAP(numberofpointdofs),stat=err)
2160  IF(err/=0) CALL flagerror("Could not allocate Neumann point DOF global to local mapping.",err,error,*999)
2161  ALLOCATE(localdofnumbers(0:rowmapping%NUMBER_OF_DOMAINS-1),stat=err)
2162  IF(err/=0) CALL flagerror("Could not allocate local Neumann DOF numbers.",err,error,*999)
2163  localdofnumbers=0
2164 
2165  IF(diagnostics2) THEN
2166  CALL write_string(diagnostic_output_type,"Local numbering",err,error,*999)
2167  END IF
2168  DO neumannidx=1,numberofpointdofs
2169  globaldof=boundaryconditionsneumann%setDofs(neumannidx)
2170  ! Get domain information from the RHS variable domain mapping, but set new local numbers.
2171  numberofdomains=rhsvariable%DOMAIN_MAPPING%GLOBAL_TO_LOCAL_MAP(globaldof)%NUMBER_OF_DOMAINS
2172  pointdofmapping%GLOBAL_TO_LOCAL_MAP(neumannidx)%NUMBER_OF_DOMAINS=numberofdomains
2173  ALLOCATE(pointdofmapping%GLOBAL_TO_LOCAL_MAP(neumannidx)%LOCAL_NUMBER(numberofdomains),stat=err)
2174  IF(err/=0) CALL flagerror("Could not allocate Neumann DOF global to local map local number.",err,error,*999)
2175  ALLOCATE(pointdofmapping%GLOBAL_TO_LOCAL_MAP(neumannidx)%DOMAIN_NUMBER(numberofdomains),stat=err)
2176  IF(err/=0) CALL flagerror("Could not allocate Neumann DOF global to local map domain number.",err,error,*999)
2177  ALLOCATE(pointdofmapping%GLOBAL_TO_LOCAL_MAP(neumannidx)%LOCAL_TYPE(numberofdomains),stat=err)
2178  IF(err/=0) CALL flagerror("Could not allocate Neumann DOF global to local map local type.",err,error,*999)
2179  IF(diagnostics2) THEN
2180  CALL write_string_value(diagnostic_output_type," Neumann point DOF index = ",neumannidx,err,error,*999)
2181  END IF
2182  DO domainidx=1,numberofdomains
2183  domainnumber=rhsvariable%DOMAIN_MAPPING%GLOBAL_TO_LOCAL_MAP(globaldof)%DOMAIN_NUMBER(domainidx)
2184  pointdofmapping%GLOBAL_TO_LOCAL_MAP(neumannidx)%DOMAIN_NUMBER(domainidx)=domainnumber
2185  pointdofmapping%GLOBAL_TO_LOCAL_MAP(neumannidx)%LOCAL_TYPE(domainidx)= &
2186  & rhsvariable%DOMAIN_MAPPING%GLOBAL_TO_LOCAL_MAP(globaldof)%LOCAL_TYPE(domainidx)
2187  IF(pointdofmapping%GLOBAL_TO_LOCAL_MAP(neumannidx)%LOCAL_TYPE(domainidx)==domain_local_internal.OR. &
2188  & pointdofmapping%GLOBAL_TO_LOCAL_MAP(neumannidx)%LOCAL_TYPE(domainidx)==domain_local_boundary) THEN
2189  localdofnumbers(domainnumber)=localdofnumbers(domainnumber)+1
2190  pointdofmapping%GLOBAL_TO_LOCAL_MAP(neumannidx)%LOCAL_NUMBER(domainidx)=localdofnumbers(domainnumber)
2191  IF(diagnostics2) THEN
2192  CALL write_string_value(diagnostic_output_type," Global rhs var DOF = ",globaldof,err,error,*999)
2193  CALL write_string_value(diagnostic_output_type," Domain number = ",domainnumber,err,error,*999)
2194  CALL write_string_value(diagnostic_output_type," Local type = ", &
2195  & pointdofmapping%GLOBAL_TO_LOCAL_MAP(neumannidx)%LOCAL_TYPE(domainidx),err,error,*999)
2196  CALL write_string_value(diagnostic_output_type," Local number = ",localdofnumbers(domainnumber),err,error,*999)
2197  END IF
2198  ENDIF
2199  END DO
2200  END DO
2201  !Local DOFs must be numbered before ghost DOFs, so loop though again, this time numbering GHOST DOFs
2202  IF(diagnostics2) THEN
2203  CALL write_string(diagnostic_output_type,"Ghost numbering",err,error,*999)
2204  END IF
2205  DO neumannidx=1,numberofpointdofs
2206  globaldof=boundaryconditionsneumann%setDofs(neumannidx)
2207  numberofdomains=rhsvariable%DOMAIN_MAPPING%GLOBAL_TO_LOCAL_MAP(globaldof)%NUMBER_OF_DOMAINS
2208  IF(diagnostics2) THEN
2209  CALL write_string_value(diagnostic_output_type," Neumann point DOF index = ",neumannidx,err,error,*999)
2210  END IF
2211  DO domainidx=1,numberofdomains
2212  IF(pointdofmapping%GLOBAL_TO_LOCAL_MAP(neumannidx)%LOCAL_TYPE(domainidx)==domain_local_ghost) THEN
2213  domainnumber=rhsvariable%DOMAIN_MAPPING%GLOBAL_TO_LOCAL_MAP(globaldof)%DOMAIN_NUMBER(domainidx)
2214  localdofnumbers(domainnumber)=localdofnumbers(domainnumber)+1
2215  pointdofmapping%GLOBAL_TO_LOCAL_MAP(neumannidx)%LOCAL_NUMBER(domainidx)=localdofnumbers(domainnumber)
2216  IF(diagnostics2) THEN
2217  CALL write_string_value(diagnostic_output_type," Global rhs var DOF = ",globaldof,err,error,*999)
2218  CALL write_string_value(diagnostic_output_type," Domain number = ",domainnumber,err,error,*999)
2219  CALL write_string_value(diagnostic_output_type," Local number = ",localdofnumbers(domainnumber),err,error,*999)
2220  END IF
2221  ENDIF
2222  END DO
2223  END DO
2224 
2225  CALL domain_mappings_local_from_global_calculate(pointdofmapping,err,error,*999)
2226 
2227  CALL distributed_matrix_create_start(rowmapping,pointdofmapping,boundaryconditionsneumann%integrationMatrix,err,error,*999)
2228  SELECT CASE(boundaryconditionsvariable%BOUNDARY_CONDITIONS%neumannMatrixSparsity)
2229  CASE(boundary_condition_sparse_matrices)
2230  ! Work out integration matrix sparsity structure
2231  ! For a single process, compressed column would be more memory efficient, but with
2232  ! multiple processes the number of Neumann point DOFs could be more than the number
2233  ! of local row DOFs, and multiplying a compressed row matrix by a vector is faster,
2234  ! so we will use compressed row storage
2235  ALLOCATE(rowindices(rowmapping%TOTAL_NUMBER_OF_LOCAL+1),stat=err)
2236  IF(err/=0) CALL flagerror("Could not allocate Neumann integration matrix column indices.",err,error,*999)
2237  ! We don't know the number of non zeros before hand, so use a list to keep track of column indices
2238  NULLIFY(columnindiceslist)
2239  CALL list_create_start(columnindiceslist,err,error,*999)
2240  CALL list_data_type_set(columnindiceslist,list_intg_type,err,error,*999)
2241  CALL list_create_finish(columnindiceslist,err,error,*999)
2242  ! Stores the column indices for the current row
2243  NULLIFY(rowcolumnindiceslist)
2244  CALL list_create_start(rowcolumnindiceslist,err,error,*999)
2245  CALL list_data_type_set(rowcolumnindiceslist,list_intg_type,err,error,*999)
2246  CALL list_mutable_set(rowcolumnindiceslist,.true.,err,error,*999)
2247  CALL list_create_finish(rowcolumnindiceslist,err,error,*999)
2248  rowindices(1)=1
2249 
2250  DO localdof=1,rhsvariable%DOMAIN_MAPPING%TOTAL_NUMBER_OF_LOCAL
2251  localdofnyy=rhsvariable%DOF_TO_PARAM_MAP%DOF_TYPE(2,localdof)
2252  componentnumber=rhsvariable%DOF_TO_PARAM_MAP%NODE_DOF2PARAM_MAP(4,localdofnyy)
2253  ! Get topology for finding faces/lines
2254  topology=>rhsvariable%COMPONENTS(componentnumber)%DOMAIN%TOPOLOGY
2255  IF(.NOT.ASSOCIATED(topology)) THEN
2256  CALL flagerror("Field component topology is not associated.",err,error,*999)
2257  END IF
2258 
2259  SELECT CASE(rhsvariable%COMPONENTS(componentnumber)%INTERPOLATION_TYPE)
2260  CASE(field_node_based_interpolation)
2261  nodenumber=rhsvariable%DOF_TO_PARAM_MAP%NODE_DOF2PARAM_MAP(3,localdofnyy)
2262  IF(.NOT.ASSOCIATED(topology%NODES%NODES)) THEN
2263  CALL flagerror("Topology nodes are not associated.",err,error,*999)
2264  END IF
2265  IF(topology%NODES%NODES(nodenumber)%BOUNDARY_NODE) THEN
2266  SELECT CASE(rhsvariable%COMPONENTS(componentnumber)%DOMAIN%NUMBER_OF_DIMENSIONS)
2267  CASE(1)
2268  ! Only one column used, as this is the same as setting an integrated
2269  ! value so no other DOFs are affected
2270  globaldof=rhsvariable%DOMAIN_MAPPING%LOCAL_TO_GLOBAL_MAP(localdof)
2271  IF(boundaryconditionsvariable%CONDITION_TYPES(globaldof)==boundary_condition_neumann_point.OR. &
2272  & boundaryconditionsvariable%CONDITION_TYPES(globaldof)==boundary_condition_neumann_point_incremented) THEN
2273  ! Find the Neumann condition number
2274  neumannconditionnumber=0
2275  DO neumannidx=1,numberofpointdofs
2276  IF(boundaryconditionsneumann%setDofs(neumannidx)==globaldof) THEN
2277  neumannconditionnumber=neumannidx
2278  END IF
2279  END DO
2280  IF(neumannconditionnumber==0) THEN
2281  CALL flagerror("Could not find matching Neuamann condition number for global DOF "// &
2282  & trim(number_to_vstring(globaldof,"*",err,error))//" with Neumann condition set.",err,error,*999)
2283  ELSE
2284  CALL list_item_add(rowcolumnindiceslist,neumannconditionnumber,err,error,*999)
2285  END IF
2286  END IF
2287  CASE(2)
2288  ! Loop over all lines for this node and find any DOFs that have a Neumann point condition set
2289  DO lineidx=1,topology%NODES%NODES(nodenumber)%NUMBER_OF_NODE_LINES
2290  IF(.NOT.ALLOCATED(topology%LINES%LINES)) THEN
2291  CALL flagerror("Topology lines have not been calculated.",err,error,*999)
2292  END IF
2293  line=>topology%LINES%LINES(topology%NODES%NODES(nodenumber)%NODE_LINES(lineidx))
2294  IF(.NOT.line%BOUNDARY_LINE) cycle
2295  DO nodeidx=1,line%BASIS%NUMBER_OF_NODES
2296  columnnodenumber=line%NODES_IN_LINE(nodeidx)
2297  DO derividx=1,line%BASIS%NUMBER_OF_DERIVATIVES(nodeidx)
2298  derivativenumber=line%DERIVATIVES_IN_LINE(1,derividx,nodeidx)
2299  versionnumber=line%DERIVATIVES_IN_LINE(2,derividx,nodeidx)
2300  columndof=rhsvariable%COMPONENTS(componentnumber)%PARAM_TO_DOF_MAP%NODE_PARAM2DOF_MAP% &
2301  & nodes(columnnodenumber)%DERIVATIVES(derivativenumber)%VERSIONS(versionnumber)
2302  globaldof=rhsvariable%DOMAIN_MAPPING%LOCAL_TO_GLOBAL_MAP(columndof)
2303  IF(boundaryconditionsvariable%CONDITION_TYPES(globaldof)==boundary_condition_neumann_point.OR. &
2304  & boundaryconditionsvariable%CONDITION_TYPES(globaldof)== &
2305  & boundary_condition_neumann_point_incremented) THEN
2306  neumannconditionnumber=0
2307  DO neumannidx=1,numberofpointdofs
2308  IF(boundaryconditionsneumann%setDofs(neumannidx)==globaldof) THEN
2309  neumannconditionnumber=neumannidx
2310  END IF
2311  END DO
2312  IF(neumannconditionnumber==0) THEN
2313  CALL flagerror("Could not find matching Neuamann condition number for global DOF "// &
2314  & trim(number_to_vstring(globaldof,"*",err,error))//" with Neumann condition set.",err,error,*999)
2315  ELSE
2316  CALL list_item_add(rowcolumnindiceslist,neumannconditionnumber,err,error,*999)
2317  END IF
2318  END IF
2319  END DO
2320  END DO
2321  END DO
2322  CASE(3)
2323  ! Loop over all faces for this node and find any DOFs that have a Neumann point condition set
2324  DO faceidx=1,topology%NODES%NODES(nodenumber)%NUMBER_OF_NODE_FACES
2325  IF(.NOT.ALLOCATED(topology%faces%faces)) THEN
2326  CALL flagerror("Topology faces have not been calculated.",err,error,*999)
2327  END IF
2328  face=>topology%FACES%FACES(topology%NODES%NODES(nodenumber)%NODE_FACES(faceidx))
2329  IF(.NOT.face%BOUNDARY_FACE) cycle
2330  DO nodeidx=1,face%BASIS%NUMBER_OF_NODES
2331  columnnodenumber=face%NODES_IN_FACE(nodeidx)
2332  DO derividx=1,face%BASIS%NUMBER_OF_DERIVATIVES(nodeidx)
2333  derivativenumber=face%DERIVATIVES_IN_FACE(1,derividx,nodeidx)
2334  versionnumber=face%DERIVATIVES_IN_FACE(2,derividx,nodeidx)
2335  columndof=rhsvariable%COMPONENTS(componentnumber)%PARAM_TO_DOF_MAP%NODE_PARAM2DOF_MAP% &
2336  & nodes(columnnodenumber)%DERIVATIVES(derivativenumber)%VERSIONS(versionnumber)
2337  globaldof=rhsvariable%DOMAIN_MAPPING%LOCAL_TO_GLOBAL_MAP(columndof)
2338  IF(boundaryconditionsvariable%CONDITION_TYPES(globaldof)==boundary_condition_neumann_point.OR. &
2339  & boundaryconditionsvariable%CONDITION_TYPES(globaldof)== &
2340  & boundary_condition_neumann_point_incremented) THEN
2341  neumannconditionnumber=0
2342  DO neumannidx=1,numberofpointdofs
2343  IF(boundaryconditionsneumann%setDofs(neumannidx)==globaldof) THEN
2344  neumannconditionnumber=neumannidx
2345  END IF
2346  END DO
2347  IF(neumannconditionnumber==0) THEN
2348  CALL flagerror("Could not find matching Neuamann condition number for global DOF "// &
2349  & trim(number_to_vstring(globaldof,"*",err,error))//" with Neumann condition set.",err,error,*999)
2350  ELSE
2351  CALL list_item_add(rowcolumnindiceslist,neumannconditionnumber,err,error,*999)
2352  END IF
2353  END IF
2354  END DO
2355  END DO
2356  END DO
2357  CASE DEFAULT
2358  CALL flagerror("The dimension is invalid for point Neumann conditions",err,error,*999)
2359  END SELECT !number of dimensions
2360  END IF
2361  CASE(field_element_based_interpolation)
2362  CALL flagerror("Not implemented.",err,error,*999)
2363  CASE(field_constant_interpolation)
2364  CALL flagerror("Not implemented.",err,error,*999)
2365  CASE(field_grid_point_based_interpolation)
2366  CALL flagerror("Not implemented.",err,error,*999)
2367  CASE(field_gauss_point_based_interpolation)
2368  CALL flagerror("Not implemented.",err,error,*999)
2369  CASE DEFAULT
2370  CALL flagerror("The interpolation type of "// &
2371  & trim(number_to_vstring(rhsvariable%COMPONENTS(componentnumber) &
2372  & %INTERPOLATION_TYPE,"*",err,error))//" is invalid for component number "// &
2373  & trim(number_to_vstring(componentnumber,"*",err,error))//".", &
2374  & err,error,*999)
2375  END SELECT
2376 
2377  !Sort and remove duplicates
2378  CALL list_remove_duplicates(rowcolumnindiceslist,err,error,*999)
2379  !Now add all column DOFs in this row that use Neumann conditions to the overall column indices
2380  CALL list_appendlist(columnindiceslist,rowcolumnindiceslist,err,error,*999)
2381  CALL list_number_of_items_get(rowcolumnindiceslist,numberrowentries,err,error,*999)
2382  rowindices(localdof+1)=rowindices(localdof)+numberrowentries
2383  CALL list_clearitems(rowcolumnindiceslist,err,error,*999)
2384  END DO !local DOFs
2385 
2386  CALL list_destroy(rowcolumnindiceslist,err,error,*999)
2387  CALL list_detach_and_destroy(columnindiceslist,numbernonzeros,columnindices,err,error,*999)
2388  IF(diagnostics1) THEN
2389  CALL write_string(diagnostic_output_type,"Neumann integration matrix sparsity",err,error,*999)
2390  CALL write_string_value(diagnostic_output_type,"Number non-zeros = ", numbernonzeros,err,error,*999)
2391  CALL write_string_value(diagnostic_output_type,"Number columns = ",numberofpointdofs,err,error,*999)
2392  CALL write_string_value(diagnostic_output_type,"Number rows = ", &
2393  & rhsvariable%DOMAIN_MAPPING%TOTAL_NUMBER_OF_LOCAL,err,error,*999)
2394  CALL write_string_vector(diagnostic_output_type,1,1,numberofpointdofs+1,6,6, &
2395  & rowindices,'(" Row indices: ",6(X,I6))', '(6X,6(X,I6))',err,error,*999)
2396  CALL write_string_vector(diagnostic_output_type,1,1,numbernonzeros,6,6, &
2397  & columnindices,'(" Column indices: ",6(X,I6))', '(6X,6(X,I6))',err,error,*999)
2398  END IF
2399 
2400  CALL distributed_matrix_storage_type_set(boundaryconditionsneumann%integrationMatrix, &
2401  & distributed_matrix_compressed_row_storage_type,err,error,*999)
2402  CALL distributed_matrix_number_non_zeros_set(boundaryconditionsneumann%integrationMatrix,numbernonzeros,err,error,*999)
2403  CALL distributed_matrix_storage_locations_set(boundaryconditionsneumann%integrationMatrix, &
2404  & rowindices,columnindices(1:numbernonzeros),err,error,*999)
2405 
2406  DEALLOCATE(localdofnumbers)
2407  DEALLOCATE(rowindices)
2408  DEALLOCATE(columnindices)
2409  CASE(boundary_condition_full_matrices)
2410  CALL distributed_matrix_storage_type_set(boundaryconditionsneumann%integrationMatrix, &
2411  & distributed_matrix_block_storage_type,err,error,*999)
2412  CASE DEFAULT
2413  CALL flagerror("The Neumann matrix sparsity type of "// &
2414  & trim(number_to_vstring(boundaryconditionsvariable%BOUNDARY_CONDITIONS%neumannMatrixSparsity,"*",err,error))// &
2415  & " is invalid.",err,error,*999)
2416  END SELECT
2417 
2418  CALL distributed_matrix_create_finish(boundaryconditionsneumann%integrationMatrix,err,error,*999)
2419 
2420  !Set up vector of Neumann point values
2421  CALL distributed_vector_create_start(pointdofmapping,boundaryconditionsneumann%pointValues,err,error,*999)
2422  CALL distributed_vector_create_finish(boundaryconditionsneumann%pointValues,err,error,*999)
2423  mycomputationalnodenumber=computational_node_number_get(err,error)
2424  !Set point values vector from boundary conditions field parameter set
2425  DO neumannidx=1,numberofpointdofs
2426  globaldof=boundaryconditionsneumann%setDofs(neumannidx)
2427  IF(rhsvariable%DOMAIN_MAPPING%GLOBAL_TO_LOCAL_MAP(globaldof)%DOMAIN_NUMBER(1)==mycomputationalnodenumber) THEN
2428  localdof=rhsvariable%DOMAIN_MAPPING%GLOBAL_TO_LOCAL_MAP(globaldof)%LOCAL_NUMBER(1)
2429  ! Set point DOF vector value
2430  localneumannconditionidx=boundaryconditionsneumann%pointDofMapping%GLOBAL_TO_LOCAL_MAP(neumannidx)%LOCAL_NUMBER(1)
2431  CALL field_parameter_set_get_local_dof(rhsvariable%FIELD,rhsvariable%VARIABLE_TYPE, &
2432  & field_boundary_conditions_set_type,localdof,pointvalue,err,error,*999)
2433  CALL distributed_vector_values_set(boundaryconditionsneumann%pointValues, &
2434  & localneumannconditionidx,pointvalue,err,error,*999)
2435  END IF
2436  END DO
2437  CALL distributed_vector_update_start(boundaryconditionsneumann%pointValues,err,error,*999)
2438  CALL distributed_vector_update_finish(boundaryconditionsneumann%pointValues,err,error,*999)
2439 
2440  ELSE
2441  CALL flagerror("The boundary condition Neumann is not associated",err,error,*998)
2442  END IF
2443  ELSE
2444  CALL flagerror("Boundary conditions variable is not associated.",err,error,*998)
2445  END IF
2446 
2447  exits("BoundaryConditions_NeumannMatricesInitialise")
2448  RETURN
2449 999 IF(ALLOCATED(rowindices)) THEN
2450  DEALLOCATE(rowindices)
2451  END IF
2452  IF(ALLOCATED(columnindices)) THEN
2453  DEALLOCATE(columnindices)
2454  END IF
2455  IF(ALLOCATED(localdofnumbers)) THEN
2456  DEALLOCATE(localdofnumbers)
2457  END IF
2458  CALL boundaryconditions_neumannmatricesfinalise(boundaryconditionsvariable,dummyerr,dummyerror,*998)
2459 998 errors("BoundaryConditions_NeumannMatricesInitialise",err,error)
2460  exits("BoundaryConditions_NeumannMatricesInitialise")
2461  RETURN 1
2462 
2463  END SUBROUTINE boundaryconditions_neumannmatricesinitialise
2464 
2465  !
2466  !================================================================================================================================
2467  !
2468 
2469  !Finalise the Neumann condition information for a boundary conditions variable
2470  SUBROUTINE boundaryconditions_neumannfinalise(boundaryConditionsVariable,err,error,*)
2471 
2472  !Argument variables
2473  TYPE(boundary_conditions_variable_type), POINTER :: boundaryConditionsVariable
2474  INTEGER(INTG), INTENT(OUT) :: err
2475  TYPE(varying_string), INTENT(OUT) :: error
2476  !Local Variables
2477  TYPE(boundaryconditionsneumanntype), POINTER :: boundaryConditionsNeumann
2478 
2479  enters("BoundaryConditions_NeumannFinalise",err,error,*999)
2480 
2481  IF(ASSOCIATED(boundaryconditionsvariable)) THEN
2482  boundaryconditionsneumann=>boundaryconditionsvariable%neumannBoundaryConditions
2483  IF(ASSOCIATED(boundaryconditionsneumann)) THEN
2484  IF(ALLOCATED(boundaryconditionsneumann%setDofs)) &
2485  & DEALLOCATE(boundaryconditionsneumann%setDofs)
2486  CALL boundaryconditions_neumannmatricesfinalise(boundaryconditionsvariable,err,error,*999)
2487  DEALLOCATE(boundaryconditionsneumann)
2488  NULLIFY(boundaryconditionsvariable%neumannBoundaryConditions)
2489  END IF
2490  ELSE
2491  CALL flagerror("Boundary conditions variable is not associated.",err,error,*999)
2492  END IF
2493 
2494  exits("BoundaryConditions_NeumannFinalise")
2495  RETURN
2496 999 errorsexits("BoundaryConditions_NeumannFinalise",err,error)
2497  RETURN 1
2498  END SUBROUTINE boundaryconditions_neumannfinalise
2499 
2500  !
2501  !================================================================================================================================
2502  !
2503 
2504  !Finalise the Neumann condition matrices for a boundary conditions variable
2505  SUBROUTINE boundaryconditions_neumannmatricesfinalise(boundaryConditionsVariable,err,error,*)
2506 
2507  !Argument variables
2508  TYPE(boundary_conditions_variable_type), POINTER :: boundaryConditionsVariable
2509  INTEGER(INTG), INTENT(OUT) :: err
2510  TYPE(varying_string), INTENT(OUT) :: error
2511  !Local Variables
2512  TYPE(boundaryconditionsneumanntype), POINTER :: boundaryConditionsNeumann
2513 
2514  enters("BoundaryConditions_NeumannMatricesFinalise",err,error,*999)
2515 
2516  IF(ASSOCIATED(boundaryconditionsvariable)) THEN
2517  boundaryconditionsneumann=>boundaryconditionsvariable%neumannBoundaryConditions
2518  IF(ASSOCIATED(boundaryconditionsneumann)) THEN
2519  IF(ASSOCIATED(boundaryconditionsneumann%integrationMatrix)) &
2520  & CALL distributed_matrix_destroy(boundaryconditionsneumann%integrationMatrix,err,error,*999)
2521  IF(ASSOCIATED(boundaryconditionsneumann%pointValues)) &
2522  & CALL distributed_vector_destroy(boundaryconditionsneumann%pointValues,err,error,*999)
2523  CALL domain_mappings_mapping_finalise(boundaryconditionsneumann%pointDofMapping,err,error,*999)
2524  END IF
2525  ELSE
2526  CALL flagerror("Boundary conditions variable is not associated.",err,error,*999)
2527  END IF
2528 
2529  exits("BoundaryConditions_NeumannMatricesFinalise")
2530  RETURN
2531 999 errorsexits("BoundaryConditions_NeumannMatricesFinalise",err,error)
2532  RETURN 1
2533  END SUBROUTINE boundaryconditions_neumannmatricesfinalise
2534 
2535  !
2536  !================================================================================================================================
2537  !
2538 
2541  SUBROUTINE boundaryconditions_neumannintegrate(rhsBoundaryConditions,err,error,*)
2542 
2543  !Argument variables
2544  TYPE(boundary_conditions_variable_type), POINTER, INTENT(IN) :: rhsBoundaryConditions
2545  INTEGER(INTG), INTENT(OUT) :: err
2546  TYPE(varying_string), INTENT(OUT) :: error
2547 
2548  !Local variables
2549  INTEGER(INTG) :: componentNumber,globalDof,localDof,neumannDofIdx,myComputationalNodeNumber
2550  INTEGER(INTG) :: numberOfNeumann,neumannLocalDof,neumannDofNyy
2551  INTEGER(INTG) :: neumannGlobalDof,neumannNodeNumber,neumannLocalNodeNumber,neumannLocalDerivNumber
2552  INTEGER(INTG) :: faceIdx,lineIdx,nodeIdx,derivIdx,gaussIdx
2553  INTEGER(INTG) :: faceNumber,lineNumber
2554  INTEGER(INTG) :: ms,os,nodeNumber,derivativeNumber,versionNumber
2555  LOGICAL :: dependentGeometry
2556  REAL(DP) :: integratedValue,phim,phio
2557  TYPE(boundaryconditionsneumanntype), POINTER :: neumannConditions
2558  TYPE(basis_type), POINTER :: basis
2559  TYPE(field_type), POINTER :: geometricField
2560  TYPE(field_variable_type), POINTER :: rhsVariable
2561  TYPE(field_interpolated_point_metrics_ptr_type), POINTER :: interpolatedPointMetrics(:)
2562  TYPE(field_interpolated_point_ptr_type), POINTER :: interpolatedPoints(:)
2563  TYPE(field_interpolation_parameters_ptr_type), POINTER :: interpolationParameters(:), scalingParameters(:)
2564  TYPE(distributed_vector_type), POINTER :: integratedValues
2565  TYPE(domain_topology_type), POINTER :: topology
2566  TYPE(domain_faces_type), POINTER :: faces
2567  TYPE(domain_lines_type), POINTER :: lines
2568  TYPE(domain_face_type), POINTER :: face
2569  TYPE(domain_line_type), POINTER :: line
2570  TYPE(decomposition_type), POINTER :: decomposition
2571  TYPE(quadrature_scheme_type), POINTER :: quadratureScheme
2572 
2573  enters("BoundaryConditions_NeumannIntegrate",err,error,*999)
2574 
2575  NULLIFY(scalingparameters)
2576  NULLIFY(interpolationparameters)
2577  NULLIFY(interpolatedpoints)
2578  NULLIFY(interpolatedpointmetrics)
2579  NULLIFY(integratedvalues)
2580 
2581  neumannconditions=>rhsboundaryconditions%neumannBoundaryConditions
2582  !Check that Neumann conditions are associated, otherwise do nothing
2583  IF(ASSOCIATED(neumannconditions)) THEN
2584  rhsvariable=>rhsboundaryconditions%VARIABLE
2585  IF(.NOT.ASSOCIATED(rhsvariable)) THEN
2586  CALL flagerror("Field variable for RHS boundary conditions is not associated.",err,error,*999)
2587  END IF
2588 
2589  CALL field_geometricgeneralfieldget(rhsvariable%field,geometricfield,dependentgeometry,err,error,*999)
2590 
2591  CALL distributed_matrix_all_values_set(neumannconditions%integrationMatrix,0.0_dp,err,error,*999)
2592 
2593  numberofneumann=rhsboundaryconditions%DOF_COUNTS(boundary_condition_neumann_point) + &
2594  & rhsboundaryconditions%DOF_COUNTS(boundary_condition_neumann_point_incremented)
2595  mycomputationalnodenumber=computational_node_number_get(err,error)
2596 
2597  ! Initialise field interpolation parameters for the geometric field, which are required for the
2598  ! face/line Jacobian and scale factors
2599  CALL field_interpolation_parameters_initialise(geometricfield,interpolationparameters,err,error,*999)
2600  CALL field_interpolation_parameters_initialise(rhsvariable%field,scalingparameters,err,error,*999)
2601  CALL field_interpolated_points_initialise(interpolationparameters,interpolatedpoints,err,error,*999)
2602  CALL field_interpolatedpointsmetricsinitialise(interpolatedpoints,interpolatedpointmetrics,err,error,*999)
2603 
2604  ! Loop over all Neumann point DOFs, finding the boundary lines or faces they are on
2605  ! and integrating over them
2606  DO neumanndofidx=1,numberofneumann
2607  neumannglobaldof=neumannconditions%setDofs(neumanndofidx)
2608  IF(rhsvariable%DOMAIN_MAPPING%GLOBAL_TO_LOCAL_MAP(neumannglobaldof)%DOMAIN_NUMBER(1)==mycomputationalnodenumber) THEN
2609  neumannlocaldof=rhsvariable%DOMAIN_MAPPING%GLOBAL_TO_LOCAL_MAP(neumannglobaldof)%LOCAL_NUMBER(1)
2610  ! Get Neumann DOF component and topology for that component
2611  neumanndofnyy=rhsvariable%DOF_TO_PARAM_MAP%DOF_TYPE(2,neumannlocaldof)
2612  componentnumber=rhsvariable%DOF_TO_PARAM_MAP%NODE_DOF2PARAM_MAP(4,neumanndofnyy)
2613  topology=>rhsvariable%COMPONENTS(componentnumber)%DOMAIN%TOPOLOGY
2614  IF(.NOT.ASSOCIATED(topology)) THEN
2615  CALL flagerror("Field component topology is not associated.",err,error,*999)
2616  END IF
2617  decomposition=>rhsvariable%COMPONENTS(componentnumber)%DOMAIN%DECOMPOSITION
2618  IF(.NOT.ASSOCIATED(decomposition)) THEN
2619  CALL flagerror("Field component decomposition is not associated.",err,error,*999)
2620  END IF
2621  SELECT CASE(rhsvariable%COMPONENTS(componentnumber)%INTERPOLATION_TYPE)
2622  CASE(field_node_based_interpolation)
2623  neumannnodenumber=rhsvariable%DOF_TO_PARAM_MAP%NODE_DOF2PARAM_MAP(3,neumanndofnyy)
2624  SELECT CASE(rhsvariable%COMPONENTS(componentnumber)%DOMAIN%NUMBER_OF_DIMENSIONS)
2625  CASE(1)
2626  CALL distributed_matrix_values_set(neumannconditions%integrationMatrix,neumannlocaldof,neumanndofidx, &
2627  & 1.0_dp,err,error,*999)
2628  CASE(2)
2629  IF(.NOT.decomposition%CALCULATE_LINES) THEN
2630  CALL flagerror("Decomposition does not have lines calculated.",err,error,*999)
2631  END IF
2632  lines=>topology%LINES
2633  IF(.NOT.ASSOCIATED(lines)) THEN
2634  CALL flagerror("Mesh topology lines is not associated.",err,error,*999)
2635  END IF
2636  linesloop: DO lineidx=1,topology%NODES%NODES(neumannnodenumber)%NUMBER_OF_NODE_LINES
2637  linenumber=topology%NODES%NODES(neumannnodenumber)%NODE_LINES(lineidx)
2638  line=>topology%lines%lines(linenumber)
2639  IF(.NOT.line%BOUNDARY_LINE) &
2640  cycle linesloop
2641  basis=>line%basis
2642  IF(.NOT.ASSOCIATED(basis)) THEN
2643  CALL flagerror("Line basis is not associated.",err,error,*999)
2644  END IF
2645  neumannlocalnodenumber=0
2646  neumannlocalderivnumber=0
2647  ! Check all nodes in line to find the local numbers for the Neumann DOF, and
2648  ! make sure we don't have an integrated_only condition set on the line
2649  DO nodeidx=1,line%BASIS%NUMBER_OF_NODES
2650  nodenumber=line%NODES_IN_LINE(nodeidx)
2651  DO derividx=1,line%BASIS%NUMBER_OF_DERIVATIVES(nodeidx)
2652  derivativenumber=line%DERIVATIVES_IN_LINE(1,derividx,nodeidx)
2653  versionnumber=line%DERIVATIVES_IN_LINE(2,derividx,nodeidx)
2654  localdof=rhsvariable%COMPONENTS(componentnumber)%PARAM_TO_DOF_MAP%NODE_PARAM2DOF_MAP% &
2655  & nodes(nodenumber)%DERIVATIVES(derivativenumber)%VERSIONS(versionnumber)
2656  globaldof=rhsvariable%DOMAIN_MAPPING%LOCAL_TO_GLOBAL_MAP(localdof)
2657  IF(globaldof==neumannglobaldof) THEN
2658  neumannlocalnodenumber=nodeidx
2659  neumannlocalderivnumber=derividx
2660  ELSE IF(rhsboundaryconditions%CONDITION_TYPES(globaldof)==boundary_condition_neumann_integrated_only) THEN
2661  cycle linesloop
2662  END IF
2663  END DO
2664  END DO
2665  IF(neumannlocalnodenumber==0) THEN
2666  CALL flagerror("Could not find local Neumann node and derivative numbers in line.",err,error,*999)
2667  END IF
2668 
2669  ! Now perform actual integration
2670  quadraturescheme=>basis%QUADRATURE%QUADRATURE_SCHEME_MAP(basis_default_quadrature_scheme)%PTR
2671  IF(.NOT.ASSOCIATED(quadraturescheme)) THEN
2672  CALL flagerror("Line basis default quadrature scheme is not associated.",err,error,*999)
2673  END IF
2674  CALL field_interpolation_parameters_line_get(field_values_set_type,linenumber, &
2675  & interpolationparameters(field_u_variable_type)%ptr,err,error,*999,field_geometric_components_type)
2676  IF(rhsvariable%FIELD%SCALINGS%SCALING_TYPE/=field_no_scaling) THEN
2677  CALL field_interpolationparametersscalefactorslineget(linenumber, &
2678  & scalingparameters(field_u_variable_type)%ptr,err,error,*999)
2679  END IF
2680 
2681  DO nodeidx=1,line%BASIS%NUMBER_OF_NODES
2682  nodenumber=line%NODES_IN_LINE(nodeidx)
2683  DO derividx=1,line%BASIS%NUMBER_OF_DERIVATIVES(nodeidx)
2684  derivativenumber=line%DERIVATIVES_IN_LINE(1,derividx,nodeidx)
2685  versionnumber=line%DERIVATIVES_IN_LINE(2,derividx,nodeidx)
2686  localdof=rhsvariable%COMPONENTS(componentnumber)%PARAM_TO_DOF_MAP%NODE_PARAM2DOF_MAP% &
2687  & nodes(nodenumber)%DERIVATIVES(derivativenumber)%VERSIONS(versionnumber)
2688 
2689  ms=basis%ELEMENT_PARAMETER_INDEX(derividx,nodeidx)
2690  os=basis%ELEMENT_PARAMETER_INDEX(neumannlocalderivnumber,neumannlocalnodenumber)
2691 
2692  integratedvalue=0.0_dp
2693  ! Loop over line gauss points, adding gauss weighted terms to the integral
2694  DO gaussidx=1,quadraturescheme%NUMBER_OF_GAUSS
2695  CALL field_interpolate_gauss(first_part_deriv,basis_default_quadrature_scheme,gaussidx, &
2696  & interpolatedpoints(field_u_variable_type)%ptr,err,error,*999,field_geometric_components_type)
2697  CALL field_interpolated_point_metrics_calculate(coordinate_jacobian_line_type, &
2698  & interpolatedpointmetrics(field_u_variable_type)%ptr,err,error,*999)
2699 
2700  !Get basis function values at guass points
2701  phim=quadraturescheme%GAUSS_BASIS_FNS(ms,no_part_deriv,gaussidx)
2702  phio=quadraturescheme%GAUSS_BASIS_FNS(os,no_part_deriv,gaussidx)
2703 
2704  !Add gauss point value to total line integral
2705  integratedvalue=integratedvalue+phim*phio* &
2706  & quadraturescheme%GAUSS_WEIGHTS(gaussidx)* &
2707  & interpolatedpointmetrics(field_u_variable_type)%ptr%jacobian
2708  END DO
2709 
2710  ! Multiply by scale factors for dependent variable
2711  IF(rhsvariable%FIELD%SCALINGS%SCALING_TYPE/=field_no_scaling) THEN
2712  integratedvalue=integratedvalue* &
2713  & scalingparameters(field_u_variable_type)%ptr%SCALE_FACTORS(ms,componentnumber)* &
2714  & scalingparameters(field_u_variable_type)%ptr%SCALE_FACTORS(os,componentnumber)
2715  END IF
2716 
2717  ! Add integral term to N matrix
2718  CALL distributed_matrix_values_add(neumannconditions%integrationMatrix,localdof,neumanndofidx, &
2719  & integratedvalue,err,error,*999)
2720  END DO
2721  END DO
2722  END DO linesloop
2723  CASE(3)
2724  IF(.NOT.decomposition%CALCULATE_FACES) THEN
2725  CALL flagerror("Decomposition does not have faces calculated.",err,error,*999)
2726  END IF
2727  faces=>topology%FACES
2728  IF(.NOT.ASSOCIATED(faces)) THEN
2729  CALL flagerror("Mesh topology faces is not associated.",err,error,*999)
2730  END IF
2731  facesloop: DO faceidx=1,topology%NODES%NODES(neumannnodenumber)%NUMBER_OF_NODE_FACES
2732  facenumber=topology%NODES%NODES(neumannnodenumber)%NODE_FACES(faceidx)
2733  face=>topology%FACES%FACES(facenumber)
2734  IF(.NOT.face%BOUNDARY_FACE) &
2735  cycle facesloop
2736  basis=>face%BASIS
2737  IF(.NOT.ASSOCIATED(basis)) THEN
2738  CALL flagerror("Line face is not associated.",err,error,*999)
2739  END IF
2740  neumannlocalnodenumber=0
2741  neumannlocalderivnumber=0
2742  ! Check all nodes in the face to find the local numbers for the Neumann DOF, and
2743  ! make sure we don't have an integrated_only condition set on the face
2744  DO nodeidx=1,basis%NUMBER_OF_NODES
2745  nodenumber=face%NODES_IN_FACE(nodeidx)
2746  DO derividx=1,basis%NUMBER_OF_DERIVATIVES(nodeidx)
2747  derivativenumber=face%DERIVATIVES_IN_FACE(1,derividx,nodeidx)
2748  versionnumber=face%DERIVATIVES_IN_FACE(2,derividx,nodeidx)
2749  localdof=rhsvariable%COMPONENTS(componentnumber)%PARAM_TO_DOF_MAP%NODE_PARAM2DOF_MAP% &
2750  & nodes(nodenumber)%DERIVATIVES(derivativenumber)%VERSIONS(versionnumber)
2751  globaldof=rhsvariable%DOMAIN_MAPPING%LOCAL_TO_GLOBAL_MAP(localdof)
2752  IF(globaldof==neumannglobaldof) THEN
2753  neumannlocalnodenumber=nodeidx
2754  neumannlocalderivnumber=derividx
2755  ELSE IF(rhsboundaryconditions%CONDITION_TYPES(globaldof)==boundary_condition_neumann_integrated_only) THEN
2756  cycle facesloop
2757  END IF
2758  END DO
2759  END DO
2760  IF(neumannlocalnodenumber==0) THEN
2761  CALL flagerror("Could not find local Neumann node and derivative numbers in line.",err,error,*999)
2762  END IF
2763 
2764  ! Now perform actual integration
2765  quadraturescheme=>basis%QUADRATURE%QUADRATURE_SCHEME_MAP(basis_default_quadrature_scheme)%PTR
2766  IF(.NOT.ASSOCIATED(quadraturescheme)) THEN
2767  CALL flagerror("Face basis default quadrature scheme is not associated.",err,error,*999)
2768  END IF
2769  CALL field_interpolation_parameters_face_get(field_values_set_type,facenumber, &
2770  & interpolationparameters(field_u_variable_type)%ptr,err,error,*999,field_geometric_components_type)
2771  IF(rhsvariable%FIELD%SCALINGS%SCALING_TYPE/=field_no_scaling) THEN
2772  CALL field_interpolationparametersscalefactorsfaceget(facenumber, &
2773  & scalingparameters(field_u_variable_type)%ptr,err,error,*999)
2774  END IF
2775 
2776  DO nodeidx=1,basis%NUMBER_OF_NODES
2777  nodenumber=face%NODES_IN_FACE(nodeidx)
2778  DO derividx=1,basis%NUMBER_OF_DERIVATIVES(nodeidx)
2779  derivativenumber=face%DERIVATIVES_IN_FACE(1,derividx,nodeidx)
2780  versionnumber=face%DERIVATIVES_IN_FACE(2,derividx,nodeidx)
2781  localdof=rhsvariable%COMPONENTS(componentnumber)%PARAM_TO_DOF_MAP%NODE_PARAM2DOF_MAP% &
2782  & nodes(nodenumber)%DERIVATIVES(derivativenumber)%VERSIONS(versionnumber)
2783 
2784  ms=basis%ELEMENT_PARAMETER_INDEX(derividx,nodeidx)
2785  os=basis%ELEMENT_PARAMETER_INDEX(neumannlocalderivnumber,neumannlocalnodenumber)
2786 
2787  integratedvalue=0.0_dp
2788  ! Loop over line gauss points, adding gauss weighted terms to the integral
2789  DO gaussidx=1,quadraturescheme%NUMBER_OF_GAUSS
2790  CALL field_interpolate_gauss(first_part_deriv,basis_default_quadrature_scheme,gaussidx, &
2791  & interpolatedpoints(field_u_variable_type)%ptr,err,error,*999,field_geometric_components_type)
2792  CALL field_interpolated_point_metrics_calculate(coordinate_jacobian_area_type, &
2793  & interpolatedpointmetrics(field_u_variable_type)%ptr,err,error,*999)
2794 
2795  !Get basis function values at guass points
2796  phim=quadraturescheme%GAUSS_BASIS_FNS(ms,no_part_deriv,gaussidx)
2797  phio=quadraturescheme%GAUSS_BASIS_FNS(os,no_part_deriv,gaussidx)
2798 
2799  !Add gauss point value to total line integral
2800  integratedvalue=integratedvalue+phim*phio* &
2801  & quadraturescheme%GAUSS_WEIGHTS(gaussidx)* &
2802  & interpolatedpointmetrics(field_u_variable_type)%ptr%jacobian
2803  END DO
2804 
2805  ! Multiply by scale factors
2806  IF(rhsvariable%FIELD%SCALINGS%SCALING_TYPE/=field_no_scaling) THEN
2807  integratedvalue=integratedvalue* &
2808  & scalingparameters(field_u_variable_type)%ptr%SCALE_FACTORS(ms,componentnumber)* &
2809  & scalingparameters(field_u_variable_type)%ptr%SCALE_FACTORS(os,componentnumber)
2810  END IF
2811 
2812  ! Add integral term to N matrix
2813  CALL distributed_matrix_values_add(neumannconditions%integrationMatrix,localdof,neumanndofidx, &
2814  & integratedvalue,err,error,*999)
2815  END DO
2816  END DO
2817  END DO facesloop
2818  CASE DEFAULT
2819  CALL flagerror("The dimension is invalid for point Neumann conditions",err,error,*999)
2820  END SELECT
2821  CASE(field_element_based_interpolation)
2822  CALL flagerror("Not implemented.",err,error,*999)
2823  CASE(field_constant_interpolation)
2824  CALL flagerror("Not implemented.",err,error,*999)
2825  CASE(field_grid_point_based_interpolation)
2826  CALL flagerror("Not implemented.",err,error,*999)
2827  CASE(field_gauss_point_based_interpolation)
2828  CALL flagerror("Not implemented.",err,error,*999)
2829  CASE DEFAULT
2830  CALL flagerror("The interpolation type of "// &
2831  & trim(number_to_vstring(rhsvariable%COMPONENTS(componentnumber) &
2832  & %INTERPOLATION_TYPE,"*",err,error))//" is invalid for component number "// &
2833  & trim(number_to_vstring(componentnumber,"*",err,error))//".", &
2834  & err,error,*999)
2835  END SELECT
2836  END IF
2837  END DO
2838 
2839  CALL distributed_matrix_update_start(neumannconditions%integrationMatrix,err,error,*999)
2840  CALL distributed_matrix_update_finish(neumannconditions%integrationMatrix,err,error,*999)
2841 
2842  CALL field_parameter_set_vector_get(rhsvariable%field,rhsvariable%variable_type,field_integrated_neumann_set_type, &
2843  & integratedvalues,err,error,*999)
2844  CALL distributed_vector_all_values_set(integratedvalues,0.0_dp,err,error,*999)
2845  ! Perform matrix multiplication, f = N q, to calculate force vector from integration matrix and point values
2846  CALL distributed_matrix_by_vector_add(distributed_matrix_vector_no_ghosts_type,1.0_dp, &
2847  & neumannconditions%integrationMatrix,neumannconditions%pointValues,integratedvalues, &
2848  & err,error,*999)
2849 
2850  CALL field_parameter_set_update_start(rhsvariable%FIELD,rhsvariable%VARIABLE_TYPE,field_integrated_neumann_set_type, &
2851  & err,error,*999)
2852  IF(diagnostics1) THEN
2853  IF(dependentgeometry) THEN
2854  CALL write_string(diagnostic_output_type," Using dependent field geometry",err,error,*999)
2855  ELSE
2856  CALL write_string(diagnostic_output_type," Using undeformed geometry",err,error,*999)
2857  END IF
2858  CALL write_string_vector(diagnostic_output_type,1,1,numberofneumann,6,6,neumannconditions%setDofs, &
2859  & '(" setDofs:",6(X,I8))', '(10X,6(X,I8))',err,error,*999)
2860  CALL write_string(diagnostic_output_type," Neumann point values",err,error,*999)
2861  CALL distributed_vector_output(diagnostic_output_type,neumannconditions%pointValues,err,error,*999)
2862  CALL write_string(diagnostic_output_type," Neumann integration matrix",err,error,*999)
2863  CALL distributed_matrix_output(diagnostic_output_type,neumannconditions%integrationMatrix,err,error,*999)
2864  CALL write_string(diagnostic_output_type," Integrated values",err,error,*999)
2865  CALL distributed_vector_output(diagnostic_output_type,integratedvalues,err,error,*999)
2866  END IF
2867  CALL field_parameter_set_update_finish(rhsvariable%FIELD,rhsvariable%VARIABLE_TYPE,field_integrated_neumann_set_type, &
2868  & err,error,*999)
2869 
2870  END IF !Neumann conditions associated
2871 
2872  exits("BoundaryConditions_NeumannIntegrate")
2873  RETURN
2874 999 errorsexits("BoundaryConditions_NeumannIntegrate",err,error)
2875  RETURN 1
2876  END SUBROUTINE boundaryconditions_neumannintegrate
2877 
2878  !
2879  !================================================================================================================================
2880  !
2881 
2883  SUBROUTINE boundaryconditions_neumannsparsitytypeset(boundaryConditions,sparsityType,err,error,*)
2884 
2885  !Argument variables
2886  INTEGER(INTG), INTENT(IN) :: sparsityType
2887  INTEGER(INTG), INTENT(OUT) :: err
2888  TYPE(varying_string), INTENT(OUT) :: error
2889  !Local Variables
2890  TYPE(boundary_conditions_type), POINTER :: boundaryConditions
2891 
2892  enters("BoundaryConditions_NeumannSparsityTypeSet",err,error,*999)
2893 
2894  IF(ASSOCIATED(boundaryconditions)) THEN
2895  SELECT CASE(sparsitytype)
2896  CASE(boundary_condition_sparse_matrices)
2897  boundaryconditions%neumannMatrixSparsity=boundary_condition_sparse_matrices
2898  CASE(boundary_condition_full_matrices)
2899  boundaryconditions%neumannMatrixSparsity=boundary_condition_full_matrices
2900  CASE DEFAULT
2901  CALL flagerror("The specified Neumann integration matrix sparsity type of "// &
2902  & trim(number_to_vstring(sparsitytype,"*",err,error))//" is invalid.",err,error,*999)
2903  END SELECT
2904  ELSE
2905  CALL flagerror("Boundary conditions are not associated.",err,error,*999)
2906  END IF
2907 
2908  exits("BoundaryConditions_NeumannSparsityTypeSet")
2909  RETURN
2910 999 errorsexits("BoundaryConditions_NeumannSparsityTypeSet",err,error)
2911  RETURN 1
2912 
2913  END SUBROUTINE boundaryconditions_neumannsparsitytypeset
2914 
2915  !
2916  !================================================================================================================================
2917  !
2918 
2920  SUBROUTINE boundary_conditions_set_node(BOUNDARY_CONDITIONS,FIELD,VARIABLE_TYPE,VERSION_NUMBER,DERIVATIVE_NUMBER, &
2921  & user_node_number,component_number,condition,VALUE,err,error,*)
2922 
2923  !Argument variables
2924  TYPE(boundary_conditions_type), POINTER :: BOUNDARY_CONDITIONS
2925  TYPE(field_type), POINTER :: FIELD
2926  INTEGER(INTG), INTENT(IN) :: VARIABLE_TYPE
2927  INTEGER(INTG), INTENT(IN) :: VERSION_NUMBER
2928  INTEGER(INTG), INTENT(IN) :: DERIVATIVE_NUMBER
2929  INTEGER(INTG), INTENT(IN) :: USER_NODE_NUMBER
2930  INTEGER(INTG), INTENT(IN) :: COMPONENT_NUMBER
2931  INTEGER(INTG), INTENT(IN) :: CONDITION
2932  REAL(DP), INTENT(IN) :: VALUE
2933  INTEGER(INTG), INTENT(OUT) :: ERR
2934  TYPE(varying_string), INTENT(OUT) :: ERROR
2935  !Local Variables
2936  INTEGER(INTG) :: local_ny,global_ny
2937  TYPE(boundary_conditions_variable_type), POINTER :: BOUNDARY_CONDITIONS_VARIABLE
2938  TYPE(field_variable_type), POINTER :: FIELD_VARIABLE
2939  TYPE(varying_string) :: LOCAL_ERROR
2940 
2941  enters("BOUNDARY_CONDITIONS_SET_NODE",err,error,*999)
2942 
2943  NULLIFY(boundary_conditions_variable)
2944  NULLIFY(field_variable)
2945 
2946  IF(ASSOCIATED(boundary_conditions)) THEN
2947  IF(boundary_conditions%BOUNDARY_CONDITIONS_FINISHED) THEN
2948  CALL flagerror("Boundary conditions have been finished.",err,error,*999)
2949  ELSE
2950  IF(ASSOCIATED(field)) THEN
2951  CALL field_component_dof_get_user_node(field,variable_type,version_number,derivative_number, &
2952  & user_node_number,component_number,local_ny,global_ny,err,error,*999)
2953  CALL field_variable_get(field,variable_type,field_variable,err,error,*999)
2954  IF(ASSOCIATED(field_variable)) THEN
2955  CALL boundary_conditions_variable_get(boundary_conditions,field_variable,boundary_conditions_variable, &
2956  & err,error,*999)
2957  IF(ASSOCIATED(boundary_conditions_variable)) THEN
2958  CALL boundaryconditions_checkinterpolationtype(condition,field,variable_type,component_number,err,error,*999)
2959  CALL boundary_conditions_set_local_dof(boundary_conditions,field,variable_type, &
2960  & local_ny,condition,VALUE,err,error,*999)
2961  ELSE
2962  local_error="The boundary conditions for variable type "//trim(number_to_vstring(variable_type,"*",err,error))// &
2963  & " has not been created."
2964  CALL flagerror(local_error,err,error,*999)
2965  ENDIF
2966  ELSE
2967  CALL flagerror("The dependent field variable is not associated",err,error,*999)
2968  ENDIF
2969  ELSE
2970  CALL flagerror("The dependent field is not associated",err,error,*999)
2971  ENDIF
2972  ENDIF
2973  ELSE
2974  CALL flagerror("Boundary conditions is not associated.",err,error,*999)
2975  ENDIF
2976 
2977  exits("BOUNDARY_CONDITIONS_SET_NODE")
2978  RETURN
2979 999 errorsexits("BOUNDARY_CONDITIONS_SET_NODE",err,error)
2980  RETURN 1
2981  END SUBROUTINE boundary_conditions_set_node
2982 
2983  !
2984  !================================================================================================================================
2985  !
2986 
2988  SUBROUTINE boundaryconditions_constraindofsequal(boundaryConditions,fieldVariable,globalDofs,coefficient,err,error,*)
2989 
2990  !Argument variables
2991  TYPE(boundary_conditions_type), POINTER, INTENT(IN) :: boundaryConditions
2992  TYPE(field_variable_type), POINTER, INTENT(IN) :: fieldVariable
2993  INTEGER(INTG), INTENT(IN) :: globalDofs(:)
2994  REAL(DP), INTENT(IN) :: coefficient
2995  INTEGER(INTG), INTENT(OUT) :: err
2996  TYPE(varying_string), INTENT(OUT) :: error
2997  !Local variables
2998  INTEGER(INTG) :: numberOfDofs,dofIdx,dofIdx2
2999 
3000  enters("BoundaryConditions_ConstrainDofsEqual",err,error,*999)
3001 
3002  numberofdofs=SIZE(globaldofs,1)
3003  IF(numberofdofs<2) THEN
3004  CALL flagerror("Cannot constrain zero or 1 DOF to be equal.",err,error,*999)
3005  END IF
3006 
3007  !Check for duplicate DOFs
3008  DO dofidx=1,numberofdofs
3009  DO dofidx2=dofidx+1,numberofdofs
3010  IF(globaldofs(dofidx)==globaldofs(dofidx2)) THEN
3011  CALL flagerror("DOF number "//trim(numbertovstring(globaldofs(dofidx),"*",err,error))// &
3012  & " is duplicated in the DOFs constrained to be equal.",err,error,*999)
3013  END IF
3014  END DO
3015  END DO
3016 
3017  !Add new DOF constraints
3018  !We set all DOFs except the first to be equal to coefficient * the first DOF
3019  !The first DOF is left unconstrained
3020  DO dofidx=2,numberofdofs
3021  CALL boundaryconditions_dofconstraintset( &
3022  & boundaryconditions,fieldvariable,globaldofs(dofidx),[globaldofs(1)],[coefficient],err,error,*999)
3023  END DO
3024 
3025  exits("BoundaryConditions_ConstrainDofsEqual")
3026  RETURN
3027 999 errorsexits("BoundaryConditions_ConstrainDofsEqual",err,error)
3028  RETURN 1
3029  END SUBROUTINE boundaryconditions_constraindofsequal
3030 
3031  !
3032  !================================================================================================================================
3033  !
3034 
3036  SUBROUTINE boundaryconditions_constrainnodedofsequal( &
3037  & boundaryconditions,field,fieldvariabletype,versionnumber,derivativenumber,component,nodes,coefficient,err,error,*)
3038 
3039  !Argument variables
3040  TYPE(boundary_conditions_type), POINTER, INTENT(IN) :: boundaryConditions
3041  TYPE(field_type), POINTER, INTENT(IN) :: field
3042  INTEGER(INTG), INTENT(IN) :: fieldVariableType
3043  INTEGER(INTG), INTENT(IN) :: versionNumber
3044  INTEGER(INTG), INTENT(IN) :: derivativeNumber
3045  INTEGER(INTG), INTENT(IN) :: component
3046  INTEGER(INTG), INTENT(IN) :: nodes(:)
3047  REAL(DP), INTENT(IN) :: coefficient
3048  INTEGER(INTG), INTENT(OUT) :: err
3049  TYPE(varying_string), INTENT(OUT) :: error
3050  !Local variables
3051  TYPE(field_variable_type), POINTER :: fieldVariable
3052  INTEGER(INTG) :: numberOfNodes, nodeIdx, dof
3053  INTEGER(INTG), ALLOCATABLE :: globalDofs(:)
3054 
3055  enters("BoundaryConditions_ConstrainNodeDofsEqual",err,error,*998)
3056 
3057  NULLIFY(fieldvariable)
3058 
3059  IF(.NOT.ASSOCIATED(boundaryconditions)) THEN
3060  CALL flagerror("Boundary conditions are not associated.",err,error,*998)
3061  END IF
3062 
3063  numberofnodes=SIZE(nodes,1)
3064  ALLOCATE(globaldofs(numberofnodes),stat=err)
3065  IF(err/=0) CALL flagerror("Could not allocate equal global DOFs array.",err,error,*998)
3066  !Get field DOFs for nodes
3067  DO nodeidx=1,numberofnodes
3068  CALL field_component_dof_get_user_node(field,fieldvariabletype,versionnumber,derivativenumber,nodes(nodeidx), &
3069  & component,dof,globaldofs(nodeidx),err,error,*999)
3070  END DO
3071  !Get the field variable and boundary conditions variable for the field
3072  CALL field_variable_get(field,fieldvariabletype,fieldvariable,err,error,*999)
3073 
3074  !Now set DOF constraint
3075  CALL boundaryconditions_constraindofsequal(boundaryconditions,fieldvariable,globaldofs,coefficient,err,error,*999)
3076 
3077  DEALLOCATE(globaldofs)
3078 
3079  exits("BoundaryConditions_ConstrainNodeDofsEqual")
3080  RETURN
3081 999 IF(ALLOCATED(globaldofs)) DEALLOCATE(globaldofs)
3082 998 errorsexits("BoundaryConditions_ConstrainNodeDofsEqual",err,error)
3083  RETURN 1
3084  END SUBROUTINE boundaryconditions_constrainnodedofsequal
3085 
3086  !
3087  !================================================================================================================================
3088  !
3089 
3091  SUBROUTINE boundaryconditions_dofconstraintset(boundaryConditions,fieldVariable,globalDof,dofs,coefficients,err,error,*)
3092 
3093  !Argument variables
3094  TYPE(boundary_conditions_type), POINTER, INTENT(IN) :: boundaryConditions
3095  TYPE(field_variable_type), POINTER, INTENT(IN) :: fieldVariable
3096  INTEGER(INTG), INTENT(IN) :: globalDof
3097  INTEGER(INTG), INTENT(IN) :: dofs(:)
3098  REAL(DP), INTENT(IN) :: coefficients(:)
3099  INTEGER(INTG), INTENT(OUT) :: err
3100  TYPE(varying_string), INTENT(OUT) :: error
3101  !Local variables
3102  INTEGER(INTG) :: numberOfDofs,dofIdx,dofIdx2
3103  TYPE(boundaryconditionsdofconstraintptrtype), ALLOCATABLE :: newConstraints(:)
3104  TYPE(boundary_conditions_variable_type), POINTER :: boundaryConditionsVariable
3105  TYPE(boundaryconditionsdofconstraintstype), POINTER :: dofConstraints
3106  TYPE(boundaryconditionsdofconstrainttype), POINTER :: dofConstraint
3107 
3108  NULLIFY(dofconstraint)
3109  NULLIFY(dofconstraints)
3110 
3111  enters("BoundaryConditions_DofConstraintSet",err,error,*998)
3112 
3113  !Check pointers for association
3114  IF(.NOT.ASSOCIATED(boundaryconditions)) THEN
3115  CALL flagerror("Boundary conditions are not associated.",err,error,*998)
3116  END IF
3117  IF(boundaryconditions%boundary_conditions_finished) THEN
3118  CALL flagerror("The boundary conditions have already been finished.",err,error,*998)
3119  END IF
3120  IF(.NOT.ASSOCIATED(fieldvariable)) THEN
3121  CALL flagerror("Field variable is not associated.",err,error,*998)
3122  END IF
3123  CALL boundary_conditions_variable_get(boundaryconditions,fieldvariable,boundaryconditionsvariable,err,error,*998)
3124  IF(.NOT.ASSOCIATED(boundaryconditionsvariable)) THEN
3125  CALL flagerror("Boundary conditions variable is not associated.",err,error,*998)
3126  END IF
3127  dofconstraints=>boundaryconditionsvariable%dofConstraints
3128  IF(.NOT.ASSOCIATED(dofconstraints)) THEN
3129  CALL flagerror("Boundary conditions DOF constraints are not associated.",err,error,*998)
3130  END IF
3131 
3132  numberofdofs=SIZE(dofs,1)
3133  IF(numberofdofs==0) THEN
3134  CALL flagerror("Empty DOFs list.",err,error,*998)
3135  ELSE IF(numberofdofs/=SIZE(coefficients,1)) THEN
3136  CALL flagerror("Length of coefficients does not match length of DOFs array.",err,error,*998)
3137  ELSE IF(numberofdofs>1) THEN
3138  CALL flagerror("Support for constraining an equations DOF to be depended on multiple "// &
3139  & "other DOFs is not yet implemented.",err,error,*998)
3140  END IF
3141 
3142  !Check for duplicate DOFs
3143  DO dofidx=1,numberofdofs
3144  DO dofidx2=dofidx+1,numberofdofs
3145  IF(dofs(dofidx)==dofs(dofidx2)) THEN
3146  CALL flagerror("DOF number "//trim(numbertovstring(dofs(dofidx),"*",err,error))// &
3147  & " is duplicated in the DOF constraint.",err,error,*998)
3148  END IF
3149  END DO
3150  END DO
3151 
3152  !Check DOFs are free
3153  DO dofidx=1,numberofdofs
3154  IF(boundaryconditionsvariable%dof_types(dofs(dofidx))/=boundary_condition_dof_free) THEN
3155  CALL flagerror("DOF number "//trim(numbertovstring(dofs(dofidx),"*",err,error))// &
3156  & " is not free in the boundary conditions.",err,error,*998)
3157  END IF
3158  END DO
3159 
3160  !Allocate new DOF constraints and copy over old constraints
3161  ALLOCATE(newconstraints(dofconstraints%numberOfConstraints+1),stat=err)
3162  IF(err/=0) CALL flagerror("Could not allocate new DOF constraints array.",err,error,*998)
3163  IF(dofconstraints%numberOfConstraints>0) THEN
3164  newconstraints(1:dofconstraints%numberOfConstraints)= &
3165  & dofconstraints%constraints(1:dofconstraints%numberOfConstraints)
3166  END IF
3167 
3168  !Set the new DOF constraint
3169  ALLOCATE(dofconstraint,stat=err)
3170  IF(err/=0) CALL flagerror("Could not allocate new DOF constraint.",err,error,*999)
3171  ALLOCATE(dofconstraint%dofs(numberofdofs),stat=err)
3172  IF(err/=0) CALL flagerror("Could not allocate constraint DOFs array.",err,error,*999)
3173  ALLOCATE(dofconstraint%coefficients(numberofdofs),stat=err)
3174  IF(err/=0) CALL flagerror("Could not allocate constraint coefficients array.",err,error,*999)
3175  dofconstraint%globalDof=globaldof
3176  dofconstraint%numberOfDofs=numberofdofs
3177  dofconstraint%dofs(1:numberofdofs)=dofs(1:numberofdofs)
3178  dofconstraint%coefficients(1:numberofdofs)=coefficients(1:numberofdofs)
3179 
3180  !Add new DOF constraint to new array
3181  newconstraints(dofconstraints%numberOfConstraints+1)%ptr=>dofconstraint
3182  !Replace old DOF constraints with new ones
3183  CALL move_alloc(newconstraints,dofconstraints%constraints)
3184  dofconstraints%numberOfConstraints=dofconstraints%numberOfConstraints+1
3185 
3186  !Set the DOF type and BC type of the constrained DOF
3187  CALL boundaryconditions_setconditiontype(boundaryconditionsvariable,globaldof,boundary_condition_linear_constraint, &
3188  & err,error,*999)
3189 
3190  exits("BoundaryConditions_DofConstraintSet")
3191  RETURN
3192 999 IF(ASSOCIATED(dofconstraint)) THEN
3193  IF(ALLOCATED(dofconstraint%dofs)) DEALLOCATE(dofconstraint%dofs)
3194  IF(ALLOCATED(dofconstraint%coefficients)) DEALLOCATE(dofconstraint%coefficients)
3195  DEALLOCATE(dofconstraint)
3196  END IF
3197  IF(ALLOCATED(newconstraints)) DEALLOCATE(newconstraints)
3198 998 errorsexits("BoundaryConditions_DofConstraintSet",err,error)
3199  RETURN 1
3200  END SUBROUTINE boundaryconditions_dofconstraintset
3201 
3202  !
3203  !================================================================================================================================
3204  !
3205 
3207  SUBROUTINE boundaryconditions_dofconstraintscreatefinish(boundaryConditionsVariable,err,error,*)
3208 
3209  !Argument variables
3210  TYPE(boundary_conditions_variable_type), POINTER :: boundaryConditionsVariable
3211  INTEGER(INTG), INTENT(OUT) :: err
3212  TYPE(varying_string), INTENT(OUT) :: error
3213  !Local variables
3214  INTEGER(INTG) :: constraintIdx,dofIdx,thisDofDomain,otherDofDomain
3215  INTEGER(INTG) :: globalDof,globalDof2,localDof,localDof2
3216  INTEGER(INTG) :: numberOfCoupledDofs
3217  INTEGER(INTG), ALLOCATABLE :: newCoupledGlobalDofs(:),newCoupledLocalDofs(:)
3218  REAL(DP), ALLOCATABLE :: newCoefficients(:)
3219  TYPE(boundaryconditionsdofconstraintstype), POINTER :: dofConstraints
3220  TYPE(boundaryconditionsdofconstrainttype), POINTER :: dofConstraint
3221  TYPE(boundaryconditionscoupleddofstype), POINTER :: dofCoupling
3222  TYPE(domain_mapping_type), POINTER :: variableDomainMapping
3223  TYPE(field_variable_type), POINTER :: fieldVariable
3224 
3225  enters("BoundaryConditions_DofConstraintsCreateFinish",err,error,*998)
3226 
3227  NULLIFY(dofcoupling)
3228 
3229  !We have a list of DOF constraints, which give the values for a field variable
3230  !DOF as a linear combination of other DOFs.
3231  !In order to be able to construct the solver matrices in the solver mapping routines,
3232  !we create a set of couplings, where a coupling is a set of field variable DOFs
3233  !mapped to a single solver row or column.
3234 
3235  IF(ASSOCIATED(boundaryconditionsvariable)) THEN
3236  fieldvariable=>boundaryconditionsvariable%variable
3237  IF(ASSOCIATED(fieldvariable)) THEN
3238  IF(ASSOCIATED(boundaryconditionsvariable%dofConstraints)) THEN
3239  dofconstraints=>boundaryconditionsvariable%dofConstraints
3240  ELSE
3241  CALL flagerror("Boundary conditions DOF constraints are not associated.",err,error,*998)
3242  END IF
3243 
3244  variabledomainmapping=>fieldvariable%domain_mapping
3245  IF(.NOT.ASSOCIATED(variabledomainmapping)) THEN
3246  CALL flagerror("Field variable domain mapping is not associated for variable type "// &
3247  & trim(numbertovstring(fieldvariable%variable_type,"*",err,error))//".",err,error,*998)
3248  END IF
3249 
3250  !Allocate an array of pointers to DOF couplings
3251  IF(dofconstraints%numberOfConstraints>0) THEN
3252  ALLOCATE(dofconstraints%dofCouplings(fieldvariable%number_of_global_dofs),stat=err)
3253  IF(err/=0) CALL flagerror( &
3254  & "Could not allocate dof constraints dof couplings array.",err,error,*998)
3255  dofconstraints%numberOfDofs=fieldvariable%number_of_global_dofs
3256  DO dofidx=1,fieldvariable%number_of_global_dofs
3257  NULLIFY(dofconstraints%dofCouplings(dofidx)%ptr)
3258  END DO
3259  END IF
3260 
3261  !Loop over all constraints
3262  DO constraintidx=1,dofconstraints%numberOfConstraints
3263  dofconstraint=>dofconstraints%constraints(constraintidx)%ptr
3264  IF(.NOT.ASSOCIATED(dofconstraint)) THEN
3265  CALL flagerror("DOF constraint number "// &
3266  & trim(numbertovstring(constraintidx,"*",err,error))// &
3267  & " is not associated.",err,error,*999)
3268  END IF
3269 
3270  globaldof=dofconstraint%globalDof
3271  localdof=variabledomainmapping%global_to_local_map(globaldof)%local_number(1)
3272  thisdofdomain=variabledomainmapping%global_to_local_map(globaldof)%domain_number(1)
3273 
3274  !Check that the constrained DOFs are still set to be constrained, as
3275  !subsequently setting a boundary condition would change the DOF type but
3276  !not update the DOF constraints structure.
3277  IF(boundaryconditionsvariable%dof_types(globaldof)/=boundary_condition_dof_constrained) THEN
3278  CALL flagerror("Global DOF number "//trim(numbertovstring(globaldof,"*",err,error))// &
3279  & " is part of a linear constraint but the DOF type has been changed"// &
3280  & " by applying a boundary condition.",err,error,*999)
3281  END IF
3282 
3283  DO dofidx=1,dofconstraint%numberOfDofs
3284  globaldof2=dofconstraint%dofs(dofidx)
3285  localdof2=variabledomainmapping%global_to_local_map(globaldof2)%local_number(1)
3286  !Check a Dirichlet conditions hasn't also been set on this DOF
3287  IF(boundaryconditionsvariable%dof_types(globaldof2)/=boundary_condition_dof_free) THEN
3288  CALL flagerror("A Dirichlet boundary condition has been set on DOF number "// &
3289  & trim(numbertovstring(globaldof2,"*",err,error))// &
3290  & " which is part of a linear constraint.",err,error,*999)
3291  END IF
3292 
3293  !Check we don't have DOF constraints that are split over domains
3294  !\todo Implement support for DOF constraints that are split over domains
3295  IF(variabledomainmapping%number_of_domains>1) THEN
3296  otherdofdomain=variabledomainmapping%global_to_local_map(globaldof2)%domain_number(1)
3297  IF(thisdofdomain/=otherdofdomain) THEN
3298  CALL flagerror("An equal DOF constraint is split over multiple domains, "// &
3299  & "support for this has not yet been implemented.",err,error,*999)
3300  END IF
3301  END IF
3302 
3303  !Add to the DOFs that are coupled with globalDof2
3304  !This might be quite inefficient if there are many dofs mapped to a single row/column
3305  !due to the reallocation at each step
3306  IF(ASSOCIATED(dofconstraints%dofCouplings(globaldof2)%ptr)) THEN
3307  numberofcoupleddofs=dofconstraints%dofCouplings(globaldof2)%ptr%numberOfDofs
3308  ALLOCATE(newcoupledglobaldofs(numberofcoupleddofs+1),stat=err)
3309  IF(err/=0) CALL flagerror("Could not allocate new DOF coupling global DOFs.",err,error,*999)
3310  ALLOCATE(newcoupledlocaldofs(numberofcoupleddofs+1),stat=err)
3311  IF(err/=0) CALL flagerror("Could not allocate new DOF coupling local DOFs.",err,error,*999)
3312  ALLOCATE(newcoefficients(numberofcoupleddofs+1),stat=err)
3313  IF(err/=0) CALL flagerror("Could not allocate new DOF coupling values.",err,error,*999)
3314  newcoupledglobaldofs(1:numberofcoupleddofs)=dofcoupling%globalDofs(1:numberofcoupleddofs)
3315  newcoupledlocaldofs(1:numberofcoupleddofs)=dofcoupling%localDofs(1:numberofcoupleddofs)
3316  newcoefficients(1:numberofcoupleddofs)=dofcoupling%coefficients(1:numberofcoupleddofs)
3317  ELSE
3318  !Set up a a new dofCoupling and set globalDof2 as the first DOF
3319  ALLOCATE(dofconstraints%dofCouplings(globaldof2)%ptr,stat=err)
3320  IF(err/=0) CALL flagerror("Could not allocate new DOF coupling type.",err,error,*999)
3321  ALLOCATE(newcoupledglobaldofs(2),stat=err)
3322  IF(err/=0) CALL flagerror("Could not allocate new DOF coupling global DOFs.",err,error,*999)
3323  ALLOCATE(newcoupledlocaldofs(2),stat=err)
3324  IF(err/=0) CALL flagerror("Could not allocate new DOF coupling local DOFs.",err,error,*999)
3325  ALLOCATE(newcoefficients(2),stat=err)
3326  IF(err/=0) CALL flagerror("Could not allocate new DOF coupling values.",err,error,*999)
3327  newcoupledglobaldofs(1)=globaldof2
3328  newcoupledlocaldofs(1)=localdof2
3329  newcoefficients(1)=1.0_dp
3330  numberofcoupleddofs=1
3331  END IF
3332  dofcoupling=>dofconstraints%dofCouplings(globaldof2)%ptr
3333  newcoupledglobaldofs(numberofcoupleddofs+1)=globaldof
3334  newcoupledlocaldofs(numberofcoupleddofs+1)=localdof
3335  newcoefficients(numberofcoupleddofs+1)=dofconstraint%coefficients(dofidx)
3336  CALL move_alloc(newcoupledglobaldofs,dofcoupling%globalDofs)
3337  CALL move_alloc(newcoupledlocaldofs,dofcoupling%localDofs)
3338  CALL move_alloc(newcoefficients,dofcoupling%coefficients)
3339  dofcoupling%numberOfDofs=numberofcoupleddofs+1
3340  END DO !dofIdx
3341  END DO !constraintIdx
3342  ELSE
3343  CALL flagerror("Field variable is not associated for this boundary conditions variable",err,error,*999)
3344  ENDIF
3345  ELSE
3346  CALL flagerror("Boundary conditions variable is not associated.",err,error,*999)
3347  END IF
3348 
3349  exits("BoundaryConditions_DofConstraintsCreateFinish")
3350  RETURN
3351 999 IF(ALLOCATED(newcoupledglobaldofs)) DEALLOCATE(newcoupledglobaldofs)
3352  IF(ALLOCATED(newcoupledlocaldofs)) DEALLOCATE(newcoupledlocaldofs)
3353  IF(ALLOCATED(newcoefficients)) DEALLOCATE(newcoefficients)
3354  CALL boundaryconditions_dofconstraintsfinalise(dofconstraints,err,error,*998)
3355 998 errors("BoundaryConditions_DofConstraintsCreateFinish",err,error)
3356  exits("BoundaryConditions_DofConstraintsCreateFinish")
3357  RETURN 1
3358 
3359  END SUBROUTINE boundaryconditions_dofconstraintscreatefinish
3360 
3361  !
3362  !================================================================================================================================
3363  !
3364 
3366  SUBROUTINE boundaryconditions_dofconstraintsfinalise(dofConstraints,err,error,*)
3367 
3368  !Argument variables
3369  TYPE(boundaryconditionsdofconstraintstype), POINTER :: dofConstraints
3370  INTEGER(INTG), INTENT(OUT) :: err
3371  TYPE(varying_string), INTENT(OUT) :: error
3372  !Local variables
3373  INTEGER(INTG) :: constraintIdx,dofIdx
3374 
3375  enters("BoundaryConditions_DofConstraintsFinalise",err,error,*999)
3376 
3377  IF(ASSOCIATED(dofconstraints)) THEN
3378  IF(ALLOCATED(dofconstraints%constraints)) THEN
3379  DO constraintidx=1,dofconstraints%numberOfConstraints
3380  IF(ASSOCIATED(dofconstraints%constraints(constraintidx)%ptr)) THEN
3381  IF(ALLOCATED(dofconstraints%constraints(constraintidx)%ptr%dofs)) THEN
3382  DEALLOCATE(dofconstraints%constraints(constraintidx)%ptr%dofs)
3383  END IF
3384  IF(ALLOCATED(dofconstraints%constraints(constraintidx)%ptr%coefficients)) THEN
3385  DEALLOCATE(dofconstraints%constraints(constraintidx)%ptr%coefficients)
3386  END IF
3387  DEALLOCATE(dofconstraints%constraints(constraintidx)%ptr)
3388  END IF
3389  END DO
3390  DEALLOCATE(dofconstraints%constraints)
3391  END IF
3392  IF(ALLOCATED(dofconstraints%dofCouplings)) THEN
3393  DO dofidx=1,dofconstraints%numberOfDofs
3394  IF(ASSOCIATED(dofconstraints%dofCouplings(dofidx)%ptr)) THEN
3395  IF(ALLOCATED(dofconstraints%dofCouplings(dofidx)%ptr%globalDofs)) THEN
3396  DEALLOCATE(dofconstraints%dofCouplings(dofidx)%ptr%globalDofs)
3397  END IF
3398  IF(ALLOCATED(dofconstraints%dofCouplings(dofidx)%ptr%localDofs)) THEN
3399  DEALLOCATE(dofconstraints%dofCouplings(dofidx)%ptr%localDofs)
3400  END IF
3401  IF(ALLOCATED(dofconstraints%dofCouplings(dofidx)%ptr%coefficients)) THEN
3402  DEALLOCATE(dofconstraints%dofCouplings(dofidx)%ptr%coefficients)
3403  END IF
3404  END IF
3405  END DO
3406  DEALLOCATE(dofconstraints%dofCouplings)
3407  END IF
3408  ELSE
3409  CALL flagerror("dofConstraints pointer is not associated.",err,error,*999)
3410  END IF
3411 
3412  exits("BoundaryConditions_DofConstraintsFinalise")
3413  RETURN
3414 999 errorsexits("BoundaryConditions_DofConstraintsFinalise",err,error)
3415  RETURN 1
3416 
3417  END SUBROUTINE boundaryconditions_dofconstraintsfinalise
3418 
3419  !
3420  !================================================================================================================================
3421  !
3422 
3424  SUBROUTINE boundaryconditions_dofconstraintsinitialise(dofConstraints,err,error,*)
3425 
3426  !Argument variables
3427  TYPE(boundaryconditionsdofconstraintstype), POINTER :: dofConstraints
3428  INTEGER(INTG), INTENT(OUT) :: err
3429  TYPE(varying_string), INTENT(OUT) :: error
3430 
3431  enters("BoundaryConditions_DofConstraintsInitialise",err,error,*999)
3432 
3433  IF(ASSOCIATED(dofconstraints)) THEN
3434  dofconstraints%numberOfConstraints=0
3435  dofconstraints%numberOfDofs=0
3436  ELSE
3437  CALL flagerror("dofConstraints pointer is not associated.",err,error,*999)
3438  END IF
3439 
3440  exits("BoundaryConditions_DofConstraintsInitialise")
3441  RETURN
3442 999 errorsexits("BoundaryConditions_DofConstraintsInitialise",err,error)
3443  RETURN 1
3444 
3445  END SUBROUTINE boundaryconditions_dofconstraintsinitialise
3446 
3447  !
3448  !================================================================================================================================
3449  !
3450 
3452  SUBROUTINE boundary_conditions_variable_finalise(BOUNDARY_CONDITIONS_VARIABLE,ERR,ERROR,*)
3453 
3454  !Argument variables
3455  TYPE(boundary_conditions_variable_type), POINTER :: BOUNDARY_CONDITIONS_VARIABLE
3456  INTEGER(INTG), INTENT(OUT) :: ERR
3457  TYPE(varying_string), INTENT(OUT) :: ERROR
3458  !Local Variables
3459  TYPE(boundary_conditions_dirichlet_type), POINTER :: BOUNDARY_CONDITIONS_DIRICHLET
3460 
3461  enters("BOUNDARY_CONDITIONS_VARIABLE_FINALISE",err,error,*999)
3462 
3463  IF(ASSOCIATED(boundary_conditions_variable)) THEN
3464  IF(ALLOCATED(boundary_conditions_variable%CONDITION_TYPES)) &
3465  & DEALLOCATE(boundary_conditions_variable%CONDITION_TYPES)
3466  IF(ALLOCATED(boundary_conditions_variable%DOF_TYPES)) &
3467  & DEALLOCATE(boundary_conditions_variable%DOF_TYPES)
3468  IF(ASSOCIATED(boundary_conditions_variable%DIRICHLET_BOUNDARY_CONDITIONS)) THEN
3469  boundary_conditions_dirichlet=>boundary_conditions_variable%DIRICHLET_BOUNDARY_CONDITIONS
3470  CALL boundaryconditions_sparsityindicesarrayfinalise(boundary_conditions_dirichlet% &
3471  & linear_sparsity_indices,err,error,*999)
3472  CALL boundaryconditions_sparsityindicesarrayfinalise(boundary_conditions_dirichlet% &
3473  & dynamic_sparsity_indices,err,error,*999)
3474  IF(ALLOCATED(boundary_conditions_dirichlet%DIRICHLET_DOF_INDICES)) THEN
3475  DEALLOCATE(boundary_conditions_dirichlet%DIRICHLET_DOF_INDICES)
3476  ENDIF
3477  DEALLOCATE(boundary_conditions_dirichlet)
3478  ENDIF
3479  CALL boundaryconditions_neumannfinalise(boundary_conditions_variable,err,error,*999)
3480  IF(ASSOCIATED(boundary_conditions_variable%PRESSURE_INCREMENTED_BOUNDARY_CONDITIONS)) &
3481  & DEALLOCATE(boundary_conditions_variable%PRESSURE_INCREMENTED_BOUNDARY_CONDITIONS)
3482  IF(ASSOCIATED(boundary_conditions_variable%dofConstraints)) THEN
3483  CALL boundaryconditions_dofconstraintsfinalise(boundary_conditions_variable%dofConstraints,err,error,*999)
3484  DEALLOCATE(boundary_conditions_variable%dofConstraints)
3485  END IF
3486  DEALLOCATE(boundary_conditions_variable)
3487  ENDIF
3488 
3489  exits("BOUNDARY_CONDITIONS_VARIABLE_FINALISE")
3490  RETURN
3491 999 errorsexits("BOUNDARY_CONDITIONS_VARIABLE_FINALISE",err,error)
3492  RETURN 1
3493  END SUBROUTINE boundary_conditions_variable_finalise
3494 
3495  !
3496  !================================================================================================================================
3497  !
3498 
3500  SUBROUTINE boundaryconditions_sparsityindicesarrayfinalise(SPARSITY_INDICES_ARRAY,ERR,ERROR,*)
3501 
3502  !Argument variables
3503  TYPE(boundary_conditions_sparsity_indices_ptr_type), ALLOCATABLE :: SPARSITY_INDICES_ARRAY(:,:)
3504  INTEGER(INTG), INTENT(OUT) :: ERR
3505  TYPE(varying_string), INTENT(OUT) :: ERROR
3506  !Local Variables
3507  INTEGER(INTG) :: equ_set_idx, equ_matrix_idx
3508  TYPE(boundary_conditions_sparsity_indices_type), POINTER :: SPARSITY_INDICES
3509 
3510  enters("BoundaryConditions_SparsityIndicesArrayFinalise",err,error,*999)
3511 
3512  IF (ALLOCATED(sparsity_indices_array)) THEN
3513  DO equ_set_idx=1,SIZE(sparsity_indices_array,1)
3514  DO equ_matrix_idx=1,SIZE(sparsity_indices_array,2)
3515  sparsity_indices=>sparsity_indices_array(equ_set_idx,equ_matrix_idx)%PTR
3516  IF(ASSOCIATED(sparsity_indices)) THEN
3517  IF(ALLOCATED(sparsity_indices%SPARSE_ROW_INDICES)) THEN
3518  DEALLOCATE(sparsity_indices%SPARSE_ROW_INDICES)
3519  ENDIF
3520  IF(ALLOCATED(sparsity_indices%SPARSE_COLUMN_INDICES)) THEN
3521  DEALLOCATE(sparsity_indices%SPARSE_COLUMN_INDICES)
3522  ENDIF
3523  DEALLOCATE(sparsity_indices)
3524  ENDIF
3525  ENDDO
3526  ENDDO
3527  DEALLOCATE(sparsity_indices_array)
3528  ENDIF
3529 
3530  exits("BoundaryConditions_SparsityIndicesArrayFinalise")
3531  RETURN
3532 999 errors("BoundaryConditions_SparsityIndicesArrayFinalise",err,error)
3533  exits("BoundaryConditions_SparsityIndicesArrayFinalise")
3534  RETURN 1
3535 
3536  END SUBROUTINE boundaryconditions_sparsityindicesarrayfinalise
3537 
3538  !
3539  !================================================================================================================================
3540  !
3541 
3543  SUBROUTINE boundary_conditions_variable_initialise(BOUNDARY_CONDITIONS,FIELD_VARIABLE,ERR,ERROR,*)
3544 
3545  !Argument variables
3546  TYPE(boundary_conditions_type), POINTER :: BOUNDARY_CONDITIONS
3547  TYPE(field_variable_type), POINTER :: FIELD_VARIABLE
3548  INTEGER(INTG), INTENT(OUT) :: ERR
3549  TYPE(varying_string), INTENT(OUT) :: ERROR
3550  !Local Variables
3551  INTEGER(INTG) :: DUMMY_ERR,variable_idx
3552  TYPE(domain_mapping_type), POINTER :: VARIABLE_DOMAIN_MAPPING
3553  TYPE(varying_string) :: DUMMY_ERROR
3554  TYPE(boundary_conditions_variable_ptr_type), ALLOCATABLE :: NEW_BOUNDARY_CONDITIONS_VARIABLES(:)
3555  TYPE(boundary_conditions_variable_type), POINTER :: BOUNDARY_CONDITIONS_VARIABLE
3556 
3557  enters("BOUNDARY_CONDITIONS_VARIABLE_INITIALISE",err,error,*998)
3558 
3559  IF(ASSOCIATED(boundary_conditions)) THEN
3560  IF(ASSOCIATED(field_variable)) THEN
3561  variable_domain_mapping=>field_variable%DOMAIN_MAPPING
3562  IF(ASSOCIATED(variable_domain_mapping)) THEN
3563  NULLIFY(boundary_conditions_variable)
3564  !Check if boundary conditions variable has already been added, if so then we don't do anything as different equations
3565  !sets can have the same dependent field variables and will both want to add the variable
3566  CALL boundary_conditions_variable_get(boundary_conditions,field_variable,boundary_conditions_variable,err,error,*999)
3567  IF(.NOT.ASSOCIATED(boundary_conditions_variable)) THEN
3568  ALLOCATE(new_boundary_conditions_variables(boundary_conditions%NUMBER_OF_BOUNDARY_CONDITIONS_VARIABLES+1),stat=err)
3569  IF(err/=0) CALL flagerror("Could not allocate new boundary conditions variables array.",err,error,*998)
3570  IF(ALLOCATED(boundary_conditions%BOUNDARY_CONDITIONS_VARIABLES)) THEN
3571  DO variable_idx=1,boundary_conditions%NUMBER_OF_BOUNDARY_CONDITIONS_VARIABLES
3572  new_boundary_conditions_variables(variable_idx)%PTR=> &
3573  & boundary_conditions%BOUNDARY_CONDITIONS_VARIABLES(variable_idx)%PTR
3574  ENDDO
3575  ENDIF
3576 
3577  ALLOCATE(new_boundary_conditions_variables(boundary_conditions%NUMBER_OF_BOUNDARY_CONDITIONS_VARIABLES+1)%PTR,stat=err)
3578  IF(err/=0) CALL flagerror("Could not allocate boundary condition variable.",err,error,*998)
3579  boundary_conditions_variable=>new_boundary_conditions_variables( &
3580  & boundary_conditions%NUMBER_OF_BOUNDARY_CONDITIONS_VARIABLES+1)%PTR
3581  boundary_conditions_variable%BOUNDARY_CONDITIONS=>boundary_conditions
3582  boundary_conditions_variable%VARIABLE_TYPE=field_variable%VARIABLE_TYPE
3583  boundary_conditions_variable%VARIABLE=>field_variable
3584  ALLOCATE(boundary_conditions_variable%CONDITION_TYPES(variable_domain_mapping%NUMBER_OF_GLOBAL),stat=err)
3585  IF(err/=0) CALL flagerror("Could not allocate global boundary condition types.",err,error,*999)
3586  ALLOCATE(boundary_conditions_variable%DOF_TYPES(variable_domain_mapping%NUMBER_OF_GLOBAL),stat=err)
3587  IF(err/=0) CALL flagerror("Could not allocate global boundary condition dof types.",err,error,*999)
3588  boundary_conditions_variable%CONDITION_TYPES=boundary_condition_free
3589  boundary_conditions_variable%DOF_TYPES=boundary_condition_dof_free
3590  ALLOCATE(boundary_conditions_variable%DOF_COUNTS(max_boundary_condition_number),stat=err)
3591  IF(err/=0) CALL flagerror("Could not allocate boundary condition DOF counts array.",err,error,*999)
3592  boundary_conditions_variable%DOF_COUNTS=0
3593  NULLIFY(boundary_conditions_variable%DIRICHLET_BOUNDARY_CONDITIONS)
3594  boundary_conditions_variable%NUMBER_OF_DIRICHLET_CONDITIONS=0
3595  NULLIFY(boundary_conditions_variable%neumannBoundaryConditions)
3596  NULLIFY(boundary_conditions_variable%PRESSURE_INCREMENTED_BOUNDARY_CONDITIONS)
3597  ALLOCATE(boundary_conditions_variable%parameterSetRequired(field_number_of_set_types),stat=err)
3598  IF(err/=0) CALL flagerror("Could not allocate boundary condition parameter set required array.",err,error,*999)
3599  boundary_conditions_variable%parameterSetRequired=.false.
3600  boundary_conditions_variable%parameterSetRequired(field_values_set_type)=.true.
3601 
3602  CALL move_alloc(new_boundary_conditions_variables,boundary_conditions%BOUNDARY_CONDITIONS_VARIABLES)
3603  boundary_conditions%NUMBER_OF_BOUNDARY_CONDITIONS_VARIABLES= &
3604  & boundary_conditions%NUMBER_OF_BOUNDARY_CONDITIONS_VARIABLES+1
3605 
3606  ALLOCATE(boundary_conditions_variable%DofConstraints,stat=err)
3607  IF(err/=0) CALL flagerror("Could not allocate boundary conditions dof constraints.",err,error,*999)
3608  CALL boundaryconditions_dofconstraintsinitialise(boundary_conditions_variable%DofConstraints,err,error,*999)
3609 
3610  END IF
3611  ELSE
3612  CALL flagerror("Field variable domain mapping is not associated.",err,error,*998)
3613  ENDIF
3614  ELSE
3615  CALL flagerror("Field variable is not associated.",err,error,*998)
3616  ENDIF
3617  ELSE
3618  CALL flagerror("Boundary conditions is not associated.",err,error,*998)
3619  ENDIF
3620 
3621  exits("BOUNDARY_CONDITIONS_VARIABLE_INITIALISE")
3622  RETURN
3623 999 CALL boundary_conditions_variable_finalise(boundary_conditions_variable,dummy_err,dummy_error,*998)
3624  DEALLOCATE(new_boundary_conditions_variables)
3625 998 errorsexits("BOUNDARY_CONDITIONS_VARIABLE_INITIALISE",err,error)
3626  RETURN 1
3627  END SUBROUTINE boundary_conditions_variable_initialise
3628 
3629  !
3630  !================================================================================================================================
3631  !
3632 
3634  SUBROUTINE boundary_conditions_variable_get(BOUNDARY_CONDITIONS,FIELD_VARIABLE,BOUNDARY_CONDITIONS_VARIABLE,ERR,ERROR,*)
3635 
3636  !Argument variables
3637  TYPE(boundary_conditions_type), POINTER :: BOUNDARY_CONDITIONS
3638  TYPE(field_variable_type), POINTER :: FIELD_VARIABLE
3639  TYPE(boundary_conditions_variable_type), POINTER, INTENT(OUT) :: BOUNDARY_CONDITIONS_VARIABLE
3640  INTEGER(INTG), INTENT(OUT) :: ERR
3641  TYPE(varying_string), INTENT(OUT) :: ERROR
3642  !Local Variables
3643  INTEGER(INTG) :: variable_idx
3644  TYPE(field_variable_type), POINTER :: VARIABLE
3645  LOGICAL :: VARIABLE_FOUND
3646 
3647  enters("BOUNDARY_CONDITIONS_VARIABLE_GET",err,error,*999)
3648 
3649  NULLIFY(boundary_conditions_variable)
3650 
3651  IF(ASSOCIATED(boundary_conditions)) THEN
3652  IF(ASSOCIATED(field_variable)) THEN
3653  IF(ALLOCATED(boundary_conditions%BOUNDARY_CONDITIONS_VARIABLES)) THEN
3654  variable_found=.false.
3655  variable_idx=1
3656  DO WHILE(variable_idx<=boundary_conditions%NUMBER_OF_BOUNDARY_CONDITIONS_VARIABLES.AND..NOT.variable_found)
3657  variable=>boundary_conditions%BOUNDARY_CONDITIONS_VARIABLES(variable_idx)%PTR%VARIABLE
3658  IF(ASSOCIATED(variable)) THEN
3659  IF(variable%VARIABLE_TYPE==field_variable%VARIABLE_TYPE.AND. &
3660  & variable%FIELD%USER_NUMBER==field_variable%FIELD%USER_NUMBER) THEN
3661  IF(ASSOCIATED(variable%FIELD%REGION)) THEN
3662  IF(variable%FIELD%REGION%USER_NUMBER==field_variable%FIELD%REGION%USER_NUMBER) THEN
3663  variable_found=.true.
3664  boundary_conditions_variable=>boundary_conditions%BOUNDARY_CONDITIONS_VARIABLES(variable_idx)%PTR
3665  ENDIF
3666  ELSEIF(ASSOCIATED(variable%FIELD%INTERFACE)) THEN
3667  IF(variable%FIELD%INTERFACE%USER_NUMBER==field_variable%FIELD%INTERFACE%USER_NUMBER) THEN
3668  variable_found=.true.
3669  boundary_conditions_variable=>boundary_conditions%BOUNDARY_CONDITIONS_VARIABLES(variable_idx)%PTR
3670  ENDIF
3671  ENDIF
3672  ENDIF
3673  variable_idx=variable_idx+1
3674  ENDIF
3675  ENDDO
3676  ENDIF
3677  ELSE
3678  CALL flagerror("Field variable is not associated.",err,error,*999)
3679  ENDIF
3680  ELSE
3681  CALL flagerror("Boundary conditions is not associated.",err,error,*999)
3682  ENDIF
3683 
3684  exits("BOUNDARY_CONDITIONS_VARIABLE_GET")
3685  RETURN
3686 999 errorsexits("BOUNDARY_CONDITIONS_VARIABLE_GET",err,error)
3687  RETURN 1
3688  END SUBROUTINE boundary_conditions_variable_get
3689 
3690  !
3691  !================================================================================================================================
3692  !
3693 
3695  SUBROUTINE boundary_conditions_dirichlet_initialise(BOUNDARY_CONDITIONS_VARIABLE,ERR,ERROR,*)
3696 
3697  !Argument variables
3698  TYPE(boundary_conditions_variable_type), POINTER :: BOUNDARY_CONDITIONS_VARIABLE
3699  INTEGER(INTG), INTENT(OUT) :: ERR
3700  TYPE(varying_string), INTENT(OUT) :: ERROR
3701  !Local Variables
3702  INTEGER(INTG) :: NUMBER_OF_DIRICHLET_CONDITIONS,NUMBER_OF_LINEAR_MATRICES,NUMBER_OF_DYNAMIC_MATRICES,matrix_idx, &
3703  & MAX_NUMBER_LINEAR_MATRICES,MAX_NUMBER_DYNAMIC_MATRICES,equations_set_idx
3704  TYPE(solver_equations_type), POINTER :: SOLVER_EQUATIONS
3705  TYPE(boundary_conditions_dirichlet_type), POINTER :: BOUNDARY_CONDITIONS_DIRICHLET
3706  TYPE(equations_set_type), POINTER :: EQUATIONS_SET
3707  TYPE(equations_type), POINTER :: EQUATIONS
3708  TYPE(equations_mapping_type), POINTER :: EQUATIONS_MAPPING
3709  TYPE(equations_mapping_linear_type), POINTER :: LINEAR_MAPPING
3710  TYPE(equations_mapping_dynamic_type), POINTER :: DYNAMIC_MAPPING
3711 
3712  enters("BOUNDARY_CONDITIONS_DIRICHLET_INITIALISE",err,error,*999)
3713 
3714  IF(ASSOCIATED(boundary_conditions_variable)) THEN
3715  IF(ASSOCIATED(boundary_conditions_variable%DIRICHLET_BOUNDARY_CONDITIONS)) THEN
3716  CALL flagerror("Dirichlet boundary conditions are already associated for this boundary conditions variable." &
3717  & ,err,error,*999)
3718  ELSE
3719  ALLOCATE(boundary_conditions_variable%DIRICHLET_BOUNDARY_CONDITIONS,stat=err)
3720  IF(err/=0) CALL flagerror("Could not allocate Dirichlet Boundary Conditions",err,error,*999)
3721  boundary_conditions_dirichlet=>boundary_conditions_variable%DIRICHLET_BOUNDARY_CONDITIONS
3722  number_of_dirichlet_conditions=boundary_conditions_variable%NUMBER_OF_DIRICHLET_CONDITIONS
3723  ALLOCATE(boundary_conditions_dirichlet%DIRICHLET_DOF_INDICES(number_of_dirichlet_conditions),stat=err)
3724  IF(err/=0) CALL flagerror("Could not allocate Dirichlet DOF indices array",err,error,*999)
3725 
3726  solver_equations=>boundary_conditions_variable%BOUNDARY_CONDITIONS%SOLVER_EQUATIONS
3727  IF(ASSOCIATED(solver_equations)) THEN
3728  max_number_linear_matrices=0
3729  max_number_dynamic_matrices=0
3730  DO equations_set_idx=1,solver_equations%SOLVER_MAPPING%NUMBER_OF_EQUATIONS_SETS
3731  equations_set=>solver_equations%SOLVER_MAPPING%EQUATIONS_SETS(equations_set_idx)%PTR
3732  IF(ASSOCIATED(equations_set)) THEN
3733  equations=>equations_set%EQUATIONS
3734  IF(ASSOCIATED(equations)) THEN
3735  equations_mapping=>equations%EQUATIONS_MAPPING
3736  IF(ASSOCIATED(equations_mapping)) THEN
3737  linear_mapping=>equations_mapping%LINEAR_MAPPING
3738  dynamic_mapping=>equations_mapping%DYNAMIC_MAPPING
3739  IF(ASSOCIATED(linear_mapping)) THEN
3740  number_of_linear_matrices=linear_mapping%NUMBER_OF_LINEAR_EQUATIONS_MATRICES
3741  IF(number_of_linear_matrices>max_number_linear_matrices) &
3742  & max_number_linear_matrices=number_of_linear_matrices
3743  ENDIF
3744  IF(ASSOCIATED(dynamic_mapping)) THEN
3745  number_of_dynamic_matrices=dynamic_mapping%NUMBER_OF_DYNAMIC_EQUATIONS_MATRICES
3746  IF(number_of_dynamic_matrices>max_number_dynamic_matrices) &
3747  & max_number_dynamic_matrices=number_of_dynamic_matrices
3748  ENDIF
3749  ELSE
3750  CALL flagerror("Equations mapping is not associated.",err,error,*999)
3751  ENDIF
3752  ELSE
3753  CALL flagerror("Equations is not associated.",err,error,*999)
3754  ENDIF
3755  ELSE
3756  CALL flagerror("Equations set is not associated.",err,error,*999)
3757  ENDIF
3758  ENDDO
3759  ALLOCATE(boundary_conditions_dirichlet%LINEAR_SPARSITY_INDICES(solver_equations%SOLVER_MAPPING%NUMBER_OF_EQUATIONS_SETS, &
3760  & max_number_linear_matrices),stat=err)
3761  IF(err/=0) CALL flagerror("Could not allocate Dirichlet linear sparsity indices array",err,error,*999)
3762  ALLOCATE(boundary_conditions_dirichlet%DYNAMIC_SPARSITY_INDICES(solver_equations%SOLVER_MAPPING%NUMBER_OF_EQUATIONS_SETS,&
3763  & max_number_dynamic_matrices),stat=err)
3764  IF(err/=0) CALL flagerror("Could not allocate Dirichlet dynamic sparsity indices array",err,error,*999)
3765  DO equations_set_idx=1,solver_equations%SOLVER_MAPPING%NUMBER_OF_EQUATIONS_SETS
3766  DO matrix_idx=1,max_number_linear_matrices
3767  NULLIFY(boundary_conditions_dirichlet%LINEAR_SPARSITY_INDICES(equations_set_idx,matrix_idx)%PTR)
3768  ENDDO
3769  DO matrix_idx=1,max_number_dynamic_matrices
3770  NULLIFY(boundary_conditions_dirichlet%DYNAMIC_SPARSITY_INDICES(equations_set_idx,matrix_idx)%PTR)
3771  ENDDO
3772  ENDDO
3773  ELSE
3774  CALL flagerror("Solver equations is not associated.",err,error,*999)
3775  ENDIF
3776  ENDIF
3777  ELSE
3778  CALL flagerror("Boundary conditions variable is not associated.",err,error,*999)
3779  ENDIF
3780 
3781  exits("BOUNDARY_CONDITIONS_DIRICHLET_INITIALISE")
3782  RETURN
3783 !!TODO \todo write BOUNDARY_CONDITIONS_DIRICHLET_FINALISE
3784 999 errorsexits("BOUNDARY_CONDITIONS_DIRICHLET_INITIALISE",err,error)
3785  RETURN 1
3786 
3787  END SUBROUTINE boundary_conditions_dirichlet_initialise
3788 
3789  !
3790  !================================================================================================================================
3791  !
3792 
3794  SUBROUTINE boundaryconditions_sparsityindicesinitialise(SPARSITY_INDICES,NUMBER_OF_DIRICHLET,ERR,ERROR,*)
3795 
3796  !Argument variables
3797  TYPE(boundary_conditions_sparsity_indices_type), POINTER :: SPARSITY_INDICES
3798  INTEGER(INTG), INTENT(IN) :: NUMBER_OF_DIRICHLET
3799  INTEGER(INTG), INTENT(OUT) :: ERR
3800  TYPE(varying_string), INTENT(OUT) :: ERROR
3801  !Local Variables
3802 
3803  enters("BoundaryConditions_SparsityIndicesInitialise",err,error,*999)
3804 
3805  IF(ASSOCIATED(sparsity_indices)) THEN
3806  CALL flagerror("Sparsity Indices are already associated.",err,error,*999)
3807  ELSE
3808  ALLOCATE(sparsity_indices,stat=err)
3809  IF(err/=0) CALL flagerror("Could not allocate sparsity indicies.",err,error,*999)
3810  ALLOCATE(sparsity_indices%SPARSE_COLUMN_INDICES(number_of_dirichlet+1),stat=err)
3811  IF(err/=0) CALL flagerror("Could not allocate sparsity column indices array",err,error,*999)
3812  ENDIF
3813 
3814  exits("BoundaryConditions_SparsityIndicesInitialise")
3815  RETURN
3816 !!TODO \todo write BOUNDARY_CONDITIONS_SPARSITY_INDICES_FINALISE
3817 999 errors("BoundaryConditions_SparsityIndicesInitialise",err,error)
3818  exits("BoundaryConditions_SparsityIndicesInitialise")
3819  RETURN 1
3820 
3821  END SUBROUTINE boundaryconditions_sparsityindicesinitialise
3822 
3823  !
3824  !================================================================================================================================
3825  !
3826 
3828  SUBROUTINE boundary_conditions_pressure_incremented_initialise(BOUNDARY_CONDITIONS_VARIABLE,ERR,ERROR,*)
3829  !Argument variables
3830  TYPE(boundary_conditions_variable_type), POINTER :: BOUNDARY_CONDITIONS_VARIABLE
3831  INTEGER(INTG), INTENT(OUT) :: ERR
3832  TYPE(varying_string), INTENT(OUT) :: ERROR
3833  !Local Variables
3834  TYPE(boundary_conditions_pressure_incremented_type), POINTER :: BOUNDARY_CONDITIONS_PRESSURE_INCREMENTED
3835  INTEGER(INTG) :: NUMBER_OF_PRESSURE_INCREMENTED_CONDITIONS
3836 
3837  enters("BOUNDARY_CONDITIONS_PRESSURE_INCREMENTED_INITIALISE",err,error,*999)
3838 
3839  IF(ASSOCIATED(boundary_conditions_variable)) THEN
3840  IF(ASSOCIATED(boundary_conditions_variable%PRESSURE_INCREMENTED_BOUNDARY_CONDITIONS)) THEN
3841  CALL flagerror("Pressure incremented boundary conditions are already associated for this boundary conditions variable." &
3842  & ,err,error,*999)
3843  ELSE
3844  ALLOCATE(boundary_conditions_variable%PRESSURE_INCREMENTED_BOUNDARY_CONDITIONS,stat=err)
3845  IF(err/=0) CALL flagerror("Could not allocate Pressure incremented Boundary Conditions",err,error,*999)
3846  boundary_conditions_pressure_incremented=>boundary_conditions_variable%PRESSURE_INCREMENTED_BOUNDARY_CONDITIONS
3847  number_of_pressure_incremented_conditions=boundary_conditions_variable%DOF_COUNTS(boundary_condition_pressure_incremented)
3848  ALLOCATE(boundary_conditions_pressure_incremented%PRESSURE_INCREMENTED_DOF_INDICES &
3849  & (number_of_pressure_incremented_conditions),stat=err)
3850  IF(err/=0) CALL flagerror("Could not allocate Pressure incremented DOF indices array",err,error,*999)
3851  ENDIF
3852  ELSE
3853  CALL flagerror("Boundary conditions variable is not associated.",err,error,*999)
3854  ENDIF
3855 
3856  exits("BOUNDARY_CONDITIONS_DIRICHLET_INITIALISE")
3857  RETURN
3858 !!TODO \todo write BOUNDARY_CONDITIONS_DIRICHLET_FINALISE
3859 999 errorsexits("BOUNDARY_CONDITIONS_DIRICHLET_INITIALISE",err,error)
3860  RETURN 1
3861 
3862  END SUBROUTINE boundary_conditions_pressure_incremented_initialise
3863 
3864  !
3865  !================================================================================================================================
3866  !
3867 
3868 END MODULE boundary_conditions_routines
boundary_conditions_routines::boundaryconditions_constraindofsequal
subroutine boundaryconditions_constraindofsequal(boundaryConditions, fieldVariable, globalDofs, coefficient, err, error,)
Constrain multiple equations dependent field DOFs to be a single solver DOF in the solver equations...
Definition: boundary_condition_routines.f90:2989
boundary_conditions_routines::boundary_conditions_add_node
subroutine, public boundary_conditions_add_node(BOUNDARY_CONDITIONS, FIELD, VARIABLE_TYPE, VERSION_NUMBER, DERIVATIVE_NUMBER, USER_NODE_NUMBER, COMPONENT_NUMBER, CONDITION, VALUE, ERR, ERROR,)
Adds to the value of the specified constant and sets this as a boundary condition on the specified us...
Definition: boundary_condition_routines.f90:1999
basis_routines
This module contains all basis function routines.
Definition: basis_routines.f90:45
types::boundary_conditions_type
Contains information on the boundary conditions for the solver equations.
Definition: types.f90:1780
base_routines::enters
subroutine, public enters(NAME, ERR, ERROR,)
Records the entry into the named procedure and initialises the error code.
Definition: base_routines.f90:363
boundary_conditions_routines::boundary_condition_moved_wall
integer(intg), parameter, public boundary_condition_moved_wall
The dof is fixed as a boundary condition.
Definition: boundary_condition_routines.f90:99
coordinate_routines
This module contains all coordinate transformation and support routines.
Definition: coordinate_routines.f90:45
types::equations_type
Contains information about the equations in an equations set.
Definition: types.f90:1735
boundary_conditions_routines::boundary_condition_moved_wall_incremented
integer(intg), parameter, public boundary_condition_moved_wall_incremented
The dof is fixed as a boundary condition, to be used with load increment loop.
Definition: boundary_condition_routines.f90:109
strings::number_to_vstring
Converts a number to its equivalent varying string representation.
Definition: strings.f90:161
boundary_conditions_routines::boundary_conditions_set_local_dofs
subroutine boundary_conditions_set_local_dofs(BOUNDARY_CONDITIONS, FIELD, VARIABLE_TYPE, DOF_INDICES, CONDITIONS, VALUES, ERR, ERROR,)
Sets a boundary condition on the specified DOFs.
Definition: boundary_condition_routines.f90:1355
boundary_conditions_routines::boundary_condition_dof_free
integer(intg), parameter, public boundary_condition_dof_free
The dof is free.
Definition: boundary_condition_routines.f90:86
boundary_conditions_routines::boundaryconditions_checkinterpolationtype
subroutine boundaryconditions_checkinterpolationtype(condition, field, variableType, componentNumber, err, error,)
Checks that the specified boundary condition is appropriate for the field variable interpolation type...
Definition: boundary_condition_routines.f90:1709
boundary_conditions_routines::boundaryconditions_dofconstraintscreatefinish
subroutine boundaryconditions_dofconstraintscreatefinish(boundaryConditionsVariable, err, error,)
Finish the creation of the dof constraints for a boundary conditions variable.
Definition: boundary_condition_routines.f90:3208
boundary_conditions_routines::boundary_conditions_finalise
subroutine boundary_conditions_finalise(BOUNDARY_CONDITIONS, ERR, ERROR,)
Finalise the boundary conditions and deallocate all memory.
Definition: boundary_condition_routines.f90:752
boundary_conditions_routines::boundary_conditions_set_local_dof1
subroutine boundary_conditions_set_local_dof1(BOUNDARY_CONDITIONS, FIELD, VARIABLE_TYPE, DOF_INDEX, CONDITION, VALUE, ERR, ERROR,)
Sets a boundary condition on the specified DOF.
Definition: boundary_condition_routines.f90:1327
boundary_conditions_routines::boundary_conditions_set_element
subroutine, public boundary_conditions_set_element(BOUNDARY_CONDITIONS, FIELD, VARIABLE_TYPE, USER_ELEMENT_NUMBER, COMPONENT_NUMBER, CONDITION, VALUE, ERR, ERROR,)
Sets a boundary condition on the specified user element.
Definition: boundary_condition_routines.f90:1934
boundary_conditions_routines::boundary_condition_linear_constraint
integer(intg), parameter boundary_condition_linear_constraint
The dof is constrained to be a linear combination of other DOFs.
Definition: boundary_condition_routines.f90:113
boundary_conditions_routines::boundary_condition_neumann_integrated
integer(intg), parameter, public boundary_condition_neumann_integrated
The dof is set to a Neumann integrated boundary condition.
Definition: boundary_condition_routines.f90:102
boundary_conditions_routines::boundaryconditions_dofconstraintsinitialise
subroutine boundaryconditions_dofconstraintsinitialise(dofConstraints, err, error,)
Initialise the DOF constraints structure.
Definition: boundary_condition_routines.f90:3425
boundary_conditions_routines::boundary_conditions_variable_initialise
subroutine boundary_conditions_variable_initialise(BOUNDARY_CONDITIONS, FIELD_VARIABLE, ERR, ERROR,)
Initialise the boundary conditions variable for a variable type if that variable has not already been...
Definition: boundary_condition_routines.f90:3544
boundary_conditions_routines::boundary_condition_dirichlet
integer(intg), parameter, public boundary_condition_dirichlet
The dof is set to a Dirichlet boundary condition.
Definition: boundary_condition_routines.f90:103
types::equations_set_type
Contains information on an equations set.
Definition: types.f90:1941
boundary_conditions_routines::boundary_condition_correction_mass_increase
integer(intg), parameter, public boundary_condition_correction_mass_increase
The dof is fixed as a boundary condition, to be used with load increment loop.
Definition: boundary_condition_routines.f90:110
strings
This module contains all string manipulation and transformation routines.
Definition: strings.f90:45
types::boundary_conditions_pressure_incremented_type
Contains information on dofs associated with pressure incremented conditions.
Definition: types.f90:1816
boundary_conditions_routines::boundary_condition_neumann_point_incremented
integer(intg), parameter, public boundary_condition_neumann_point_incremented
A Neumann point boundary condition that is incremented inside a load increment control loop...
Definition: boundary_condition_routines.f90:114
boundary_conditions_routines::boundary_condition_pressure_incremented
integer(intg), parameter, public boundary_condition_pressure_incremented
The dof is a surface pressure boundary condition, to be used with load increment loop.
Definition: boundary_condition_routines.f90:108
boundary_conditions_routines::boundary_conditions_set_constant
subroutine, public boundary_conditions_set_constant(BOUNDARY_CONDITIONS, FIELD, VARIABLE_TYPE, COMPONENT_NUMBER, CONDITION, VALUE, ERR, ERROR,)
Sets a boundary condition on the specified constant.
Definition: boundary_condition_routines.f90:1076
timer
This module contains routines for timing the program.
Definition: timer_f.f90:45
boundary_conditions_routines::boundary_conditions_destroy
subroutine, public boundary_conditions_destroy(BOUNDARY_CONDITIONS, ERR, ERROR,)
Destroys boundary conditions.
Definition: boundary_condition_routines.f90:724
boundary_conditions_routines::boundaryconditions_neumanninitialise
subroutine boundaryconditions_neumanninitialise(boundaryConditionsVariable, err, error,)
Initialise the Neumann boundary conditions information.
Definition: boundary_condition_routines.f90:2065
boundary_conditions_routines::boundary_conditions_pressure_incremented_initialise
subroutine boundary_conditions_pressure_incremented_initialise(BOUNDARY_CONDITIONS_VARIABLE, ERR, ERROR,)
Initialises the pressure incremented boundary condition.
Definition: boundary_condition_routines.f90:3829
boundary_conditions_routines::boundaryconditions_sparsityindicesinitialise
subroutine boundaryconditions_sparsityindicesinitialise(SPARSITY_INDICES, NUMBER_OF_DIRICHLET, ERR, ERROR,)
Initialise Sparsity Indices type.
Definition: boundary_condition_routines.f90:3795
boundary_conditions_routines::boundaryconditions_checkequations
subroutine boundaryconditions_checkequations(boundaryConditionsVariable, err, error,)
Checks that the applied boundary conditions are supported by the equations sets in the solver equatio...
Definition: boundary_condition_routines.f90:1792
iso_varying_string
This module provides an iso_varying_string module, conformant to the API specified in ISO/IEC 1539-2:...
Definition: iso_varying_string.f90:62
linkedlist_routines
Only for integer data type for now.
Definition: linkedlist_routines.f90:2
boundary_conditions_routines::boundaryconditions_neumannintegrate
subroutine, public boundaryconditions_neumannintegrate(rhsBoundaryConditions, err, error,)
Calculates integrated Neumann condition values from point values for a boundary conditions variable a...
Definition: boundary_condition_routines.f90:2542
boundary_conditions_routines::boundaryconditions_neumannmatricesfinalise
subroutine boundaryconditions_neumannmatricesfinalise(boundaryConditionsVariable, err, error,)
Definition: boundary_condition_routines.f90:2506
boundary_conditions_routines::boundary_condition_fixed_stree
integer(intg), parameter, public boundary_condition_fixed_stree
The dof is fixed and set to values specified based on the transmission line theory at the dof...
Definition: boundary_condition_routines.f90:118
iso_varying_string::varying_string
Definition: iso_varying_string.f90:75
constants
This module contains all program wide constants.
Definition: constants.f90:45
boundary_conditions_routines::boundary_condition_fixed_inlet
integer(intg), parameter, public boundary_condition_fixed_inlet
The dof is fixed as a boundary condition.
Definition: boundary_condition_routines.f90:96
boundary_conditions_routines::max_boundary_condition_number
integer(intg), parameter max_boundary_condition_number
Definition: boundary_condition_routines.f90:121
boundary_conditions_routines::boundary_condition_fixed_nonreflecting
integer(intg), parameter, public boundary_condition_fixed_nonreflecting
The dof is fixed and set to a non-reflecting type for 1D wave propagation problems.
Definition: boundary_condition_routines.f90:116
boundary_conditions_routines::boundary_conditions_add_element
subroutine, public boundary_conditions_add_element(BOUNDARY_CONDITIONS, FIELD, VARIABLE_TYPE, USER_ELEMENT_NUMBER, COMPONENT_NUMBER, CONDITION, VALUE, ERR, ERROR,)
Adds to the value of the specified constant and sets this as a boundary condition on the specified us...
Definition: boundary_condition_routines.f90:1645
boundary_conditions_routines::boundary_condition_free
integer(intg), parameter, public boundary_condition_free
The dof is free.
Definition: boundary_condition_routines.f90:94
types::boundary_conditions_variable_type
Contains information on the boundary conditions for a dependent field variable.
Definition: types.f90:1759
boundary_conditions_routines::boundary_condition_neumann_point
integer(intg), parameter, public boundary_condition_neumann_point
The dof is set to a Neumann point boundary condition.
Definition: boundary_condition_routines.f90:101
types
This module contains all type definitions in order to avoid cyclic module references.
Definition: types.f90:70
boundary_conditions_routines::boundary_condition_fixed_outlet
integer(intg), parameter, public boundary_condition_fixed_outlet
The dof is fixed as a boundary condition.
Definition: boundary_condition_routines.f90:97
types::equations_matrices_type
Contains information on the equations matrices and vectors.
Definition: types.f90:1520
boundary_conditions_routines::boundary_conditions_add_local_dof
Adds to the value of the specified local DOF and sets this as a boundary condition on the specified l...
Definition: boundary_condition_routines.f90:137
boundary_conditions_routines::boundary_condition_fixed_wall
integer(intg), parameter, public boundary_condition_fixed_wall
The dof is fixed as a boundary condition.
Definition: boundary_condition_routines.f90:98
base_routines
This module contains all the low-level base routines e.g., all debug, control, and low-level communic...
Definition: base_routines.f90:45
types::equations_matrices_linear_type
Contains information of the linear matrices for equations matrices.
Definition: types.f90:1479
types::boundary_conditions_dirichlet_type
Contains information on dofs with associated dirichlet conditions and corresponding non-zero elements...
Definition: types.f90:1794
boundary_conditions_routines::boundary_conditions_create_finish
subroutine, public boundary_conditions_create_finish(BOUNDARY_CONDITIONS, ERR, ERROR,)
Finish the creation of boundary conditions.
Definition: boundary_condition_routines.f90:179
iso_varying_string::trim
Definition: iso_varying_string.f90:204
boundary_conditions_routines::boundary_conditions_set_node
subroutine, public boundary_conditions_set_node(BOUNDARY_CONDITIONS, FIELD, VARIABLE_TYPE, VERSION_NUMBER, DERIVATIVE_NUMBER, USER_NODE_NUMBER, COMPONENT_NUMBER, CONDITION, VALUE, ERR, ERROR,)
Sets a boundary condition on the specified user node.
Definition: boundary_condition_routines.f90:2922
boundary_conditions_routines::boundary_conditions_initialise
subroutine boundary_conditions_initialise(SOLVER_EQUATIONS, ERR, ERROR,)
Initialises the boundary conditions for an equations set.
Definition: boundary_condition_routines.f90:793
boundary_conditions_routines::boundary_condition_neumann_integrated_only
integer(intg), parameter, public boundary_condition_neumann_integrated_only
A Neumann integrated boundary condition, and no point values will be integrated over a face or line t...
Definition: boundary_condition_routines.f90:112
types::list_type
Contains information on a list.
Definition: types.f90:113
boundary_conditions_routines::boundaryconditions_neumannsparsitytypeset
subroutine, public boundaryconditions_neumannsparsitytypeset(boundaryConditions, sparsityType, err, error,)
Sets/changes the sparsity type for the Neumann integration matrices.
Definition: boundary_condition_routines.f90:2884
comp_environment
This module contains all computational environment variables.
Definition: computational_environment.f90:46
cmiss_mpi
This module contains CMISS MPI routines.
Definition: cmiss_mpi.f90:45
boundary_conditions_routines::boundary_condition_impermeable_wall
integer(intg), parameter, public boundary_condition_impermeable_wall
The dof is set such that (via penalty formulation): velocity * normal = 0.
Definition: boundary_condition_routines.f90:111
domain_mappings
This module handles all domain mappings routines.
Definition: domain_mappings.f90:45
boundary_conditions_routines::boundary_condition_fixed_cellml
integer(intg), parameter, public boundary_condition_fixed_cellml
The dof is fixed and set to values specified based on the coupled CellML solution at the dof...
Definition: boundary_condition_routines.f90:117
types::solver_equations_type
Contains information about the solver equations for a solver.
Definition: types.f90:2452
boundary_conditions_routines::boundaryconditions_dofconstraintset
subroutine boundaryconditions_dofconstraintset(boundaryConditions, fieldVariable, globalDof, dofs, coefficients, err, error,)
Constrain a DOF to be a linear combination of other DOFs.
Definition: boundary_condition_routines.f90:3092
comp_environment::computational_environment
type(computational_environment_type), target, public computational_environment
The computational environment the program is running in.
Definition: computational_environment.f90:125
boundary_conditions_routines::boundaryconditions_dofconstraintsfinalise
subroutine boundaryconditions_dofconstraintsfinalise(dofConstraints, err, error,)
Finalise the DOF constraints structure.
Definition: boundary_condition_routines.f90:3367
boundary_conditions_routines::boundary_conditions_dirichlet_initialise
subroutine boundary_conditions_dirichlet_initialise(BOUNDARY_CONDITIONS_VARIABLE, ERR, ERROR,)
Initialise dirichlet boundary conditions for a boundary conditions.
Definition: boundary_condition_routines.f90:3696
boundary_conditions_routines::boundaryconditions_neumannmatricesinitialise
subroutine boundaryconditions_neumannmatricesinitialise(boundaryConditionsVariable, err, error,)
Initialise the Neumann boundary conditions matrices and vectors. This must be done after we know whic...
Definition: boundary_condition_routines.f90:2115
boundary_conditions_routines::boundary_condition_dof_mixed
integer(intg), parameter, public boundary_condition_dof_mixed
The dof is set as a mixed boundary condition.
Definition: boundary_condition_routines.f90:88
boundary_conditions_routines::boundary_conditions_variable_finalise
subroutine boundary_conditions_variable_finalise(BOUNDARY_CONDITIONS_VARIABLE, ERR, ERROR,)
Finalise the boundary conditions variable and deallocate all memory.
Definition: boundary_condition_routines.f90:3453
distributed_matrix_vector
This module handles all distributed matrix vector routines.
Definition: distributed_matrix_vector.f90:45
interface_conditions_constants
This module defines all constants shared across interface condition routines.
Definition: interface_conditions_constants.f90:45
boundary_conditions_routines::boundary_condition_sparse_matrices
integer(intg), parameter, public boundary_condition_sparse_matrices
The matrices are stored as sparse matrices.
Definition: boundary_condition_routines.f90:126
boundary_conditions_routines
This module handles all boundary conditions routines.
Definition: boundary_condition_routines.f90:45
types::equations_matrix_type
Contains information about an equations matrix.
Definition: types.f90:1429
boundary_conditions_routines::boundaryconditions_constrainnodedofsequal
subroutine, public boundaryconditions_constrainnodedofsequal(boundaryConditions, field, fieldVariableType, versionNumber, derivativeNumber, component, nodes, coefficient, err, error,)
Constrain multiple nodal equations dependent field DOFs to be a single solver DOF in the solver equat...
Definition: boundary_condition_routines.f90:3038
boundary_conditions_routines::boundary_condition_free_wall
integer(intg), parameter, public boundary_condition_free_wall
The dof is fixed as a boundary condition.
Definition: boundary_condition_routines.f90:100
boundary_conditions_routines::boundary_condition_full_matrices
integer(intg), parameter, public boundary_condition_full_matrices
The matrices are stored as full matrices.
Definition: boundary_condition_routines.f90:127
boundary_conditions_routines::boundary_conditions_set_local_dof
Sets a boundary condition on the specified local DOF.
Definition: boundary_condition_routines.f90:143
types::field_variable_type
Contains information for a field variable defined on a field.
Definition: types.f90:1289
boundary_conditions_routines::boundaryconditions_neumannfinalise
subroutine boundaryconditions_neumannfinalise(boundaryConditionsVariable, err, error,)
Definition: boundary_condition_routines.f90:2471
boundary_conditions_routines::boundary_conditions_add_constant
subroutine, public boundary_conditions_add_constant(BOUNDARY_CONDITIONS, FIELD, VARIABLE_TYPE, COMPONENT_NUMBER, CONDITION, VALUE, ERR, ERROR,)
Adds to the value of the specified constant and sets this as a boundary condition on the specified co...
Definition: boundary_condition_routines.f90:1015
boundary_conditions_routines::boundary_condition_dof_constrained
integer(intg), parameter, public boundary_condition_dof_constrained
The dof is constrained to be a linear combination of other DOFs.
Definition: boundary_condition_routines.f90:89
types::domain_mapping_type
Contains information on the domain mappings (i.e., local and global numberings).
Definition: types.f90:904
boundary_conditions_routines::boundary_condition_fixed_fitted
integer(intg), parameter, public boundary_condition_fixed_fitted
The dof is fixed as a boundary condition to be updated from fitting data.
Definition: boundary_condition_routines.f90:115
boundary_conditions_routines::boundary_condition_fixed
integer(intg), parameter, public boundary_condition_fixed
The dof is fixed as a boundary condition.
Definition: boundary_condition_routines.f90:95
boundary_conditions_routines::boundary_conditions_add_local_dofs
subroutine boundary_conditions_add_local_dofs(BOUNDARY_CONDITIONS, FIELD, VARIABLE_TYPE, DOF_INDICES, CONDITIONS, VALUES, ERR, ERROR,)
Adds to the value of the specified DOF and sets this as a boundary condition on the specified DOFs...
Definition: boundary_condition_routines.f90:1162
equations_set_constants
This module defines all constants shared across equations set routines.
Definition: equations_set_constants.f90:45
boundary_conditions_routines::boundaryconditions_setconditiontype
subroutine boundaryconditions_setconditiontype(boundaryConditionsVariable, globalDof, condition, err, error,)
Checks the boundary condition type and sets the boundary condition type and dof type for the boundary...
Definition: boundary_condition_routines.f90:1510
boundary_conditions_routines::boundary_conditions_add_local_dof1
subroutine boundary_conditions_add_local_dof1(BOUNDARY_CONDITIONS, FIELD, VARIABLE_TYPE, DOF_INDEX, CONDITION, VALUE, ERR, ERROR,)
Adds to the value of the specified DOF and sets this as a boundary condition on the specified DOF...
Definition: boundary_condition_routines.f90:1133
boundary_conditions_routines::boundary_condition_fixed_incremented
integer(intg), parameter, public boundary_condition_fixed_incremented
The dof is a fixed boundary condition, to be used with load increment loop.
Definition: boundary_condition_routines.f90:106
lists
Implements lists of base types.
Definition: lists.f90:46
boundary_conditions_routines::boundary_conditions_variable_get
subroutine, public boundary_conditions_variable_get(BOUNDARY_CONDITIONS, FIELD_VARIABLE, BOUNDARY_CONDITIONS_VARIABLE, ERR, ERROR,)
Find the boundary conditions variable for a given field variable.
Definition: boundary_condition_routines.f90:3635
types::boundary_conditions_sparsity_indices_type
Contains information on indices of non-zero elements with associated dirichlet conditions Indices sto...
Definition: types.f90:1802
boundary_conditions_routines::boundary_condition_cauchy
integer(intg), parameter, public boundary_condition_cauchy
The dof is set to a Cauchy boundary condition.
Definition: boundary_condition_routines.f90:104
boundary_conditions_routines::boundary_condition_pressure
integer(intg), parameter, public boundary_condition_pressure
The dof is a surface pressure boundary condition.
Definition: boundary_condition_routines.f90:107
base_routines::flagerror
Flags an error condition.
Definition: base_routines.f90:240
boundary_conditions_routines::boundaryconditions_sparsityindicesarrayfinalise
subroutine boundaryconditions_sparsityindicesarrayfinalise(SPARSITY_INDICES_ARRAY, ERR, ERROR,)
Finalise an array of sparcity indices and deallocate all memory.
Definition: boundary_condition_routines.f90:3501
linkedlist_routines::linkedlist
Definition: linkedlist_routines.f90:20
boundary_conditions_routines::boundary_condition_dof_fixed
integer(intg), parameter, public boundary_condition_dof_fixed
The dof is fixed as a boundary condition.
Definition: boundary_condition_routines.f90:87
boundary_conditions_routines::boundary_conditions_create_start
subroutine, public boundary_conditions_create_start(SOLVER_EQUATIONS, BOUNDARY_CONDITIONS, ERR, ERROR,)
Start the creation of boundary conditions for the equation set.
Definition: boundary_condition_routines.f90:678
kinds
This module contains all kind definitions.
Definition: kinds.f90:45
boundary_conditions_routines::boundary_condition_robin
integer(intg), parameter, public boundary_condition_robin
The dof is set to a Robin boundary condition.
Definition: boundary_condition_routines.f90:105
types::equations_matrices_dynamic_type
Contains information of the dynamic matrices for equations matrices.
Definition: types.f90:1471
cmiss_mpi::mpi_error_check
subroutine, public mpi_error_check(ROUTINE, MPI_ERR_CODE, ERR, ERROR,)
Checks to see if an MPI error has occured during an MPI call and flags a CMISS error it if it has...
Definition: cmiss_mpi.f90:84
input_output
This module handles all formating and input and output.
Definition: input_output.f90:45