122 INTEGER(INTG),
INTENT(IN) :: SOLUTION_METHOD
123 INTEGER(INTG),
INTENT(OUT) :: ERR
128 enters(
"Darcy_EquationsSetSolutionMethodSet",err,error,*999)
130 IF(
ASSOCIATED(equations_set))
THEN 131 IF(.NOT.
ALLOCATED(equations_set%SPECIFICATION))
THEN 132 CALL flagerror(
"Equations set specification is not allocated.",err,error,*999)
133 ELSE IF(
SIZE(equations_set%SPECIFICATION,1)/=3)
THEN 134 CALL flagerror(
"Equations set specification must have three entries for a Darcy type equations set.", &
137 SELECT CASE(equations_set%SPECIFICATION(3))
143 SELECT CASE(solution_method)
147 CALL flagerror(
"Not implemented.",err,error,*999)
149 CALL flagerror(
"Not implemented.",err,error,*999)
151 CALL flagerror(
"Not implemented.",err,error,*999)
153 CALL flagerror(
"Not implemented.",err,error,*999)
155 CALL flagerror(
"Not implemented.",err,error,*999)
157 local_error=
"The specified solution method of "//
trim(
number_to_vstring(solution_method,
"*",err,error))//
" is invalid." 158 CALL flagerror(local_error,err,error,*999)
161 local_error=
"Equations set subtype of "//
trim(
number_to_vstring(equations_set%SPECIFICATION(3),
"*",err,error))// &
162 &
" is not valid for a Darcy equation type of a fluid mechanics equations set class." 163 CALL flagerror(local_error,err,error,*999)
166 CALL flagerror(
"Equations set is not associated.",err,error,*999)
169 exits(
"Darcy_EquationsSetSolutionMethodSet")
171 999 errorsexits(
"Darcy_EquationsSetSolutionMethodSet",err,error)
185 INTEGER(INTG),
INTENT(OUT) :: ERR
188 INTEGER(INTG) :: GEOMETRIC_SCALING_TYPE, GEOMETRIC_MESH_COMPONENT,NUMBER_OF_COMPONENTS,NUMBER_OF_DARCY_COMPONENTS
189 INTEGER(INTG) :: imy_matrix,Ncompartments
196 TYPE(
field_type),
POINTER :: EQUATIONS_SET_FIELD_FIELD
200 INTEGER:: DEPENDENT_FIELD_NUMBER_OF_VARIABLES, DEPENDENT_FIELD_NUMBER_OF_COMPONENTS
201 INTEGER:: DEPENDENT_FIELD_ELASTICITY_NUMBER_OF_COMPONENTS, DEPENDENT_FIELD_DARCY_NUMBER_OF_COMPONENTS
202 INTEGER:: INDEPENDENT_FIELD_NUMBER_OF_VARIABLES, INDEPENDENT_FIELD_NUMBER_OF_COMPONENTS
203 INTEGER:: NUMBER_OF_DIMENSIONS, GEOMETRIC_COMPONENT_NUMBER
204 INTEGER:: MATERIAL_FIELD_NUMBER_OF_VARIABLES, MATERIAL_FIELD_NUMBER_OF_COMPONENTS
205 INTEGER:: MESH_COMPONENT,MATERIAL_FIELD_NUMBER_OF_U_VAR_COMPONENTS,MATERIAL_FIELD_NUMBER_OF_V_VAR_COMPONENTS, &
206 & MATERIAL_FIELD_NUMBER_OF_U1_VAR_COMPONENTS
207 INTEGER:: i,component_idx
209 INTEGER(INTG) :: num_var,num_var_count
210 INTEGER(INTG) :: EQUATIONS_SET_FIELD_NUMBER_OF_VARIABLES,EQUATIONS_SET_FIELD_NUMBER_OF_COMPONENTS,NUMBER_OF_SOURCE_COMPONENTS
211 INTEGER(INTG),
POINTER :: EQUATIONS_SET_FIELD_DATA(:)
212 INTEGER(INTG),
ALLOCATABLE :: VARIABLE_TYPES(:),VARIABLE_U_TYPES(:),COUPLING_MATRIX_STORAGE_TYPE(:), &
213 & COUPLING_MATRIX_STRUCTURE_TYPE(:)
215 enters(
"DARCY_EQUATION_EQUATIONS_SET_SETUP",err,error,*999)
218 NULLIFY(equations_mapping)
219 NULLIFY(equations_matrices)
220 NULLIFY(geometric_decomposition)
221 NULLIFY(equations_equations_set_field)
222 NULLIFY(equations_set_field_field)
223 NULLIFY(equations_set_field_data)
225 IF(
ASSOCIATED(equations_set))
THEN 226 IF(.NOT.
ALLOCATED(equations_set%SPECIFICATION))
THEN 227 CALL flagerror(
"Equations set specification is not allocated.",err,error,*999)
228 ELSE IF(
SIZE(equations_set%SPECIFICATION,1)/=3)
THEN 229 CALL flagerror(
"Equations set specification must have three entries for a Darcy type equations set.", &
232 SELECT CASE(equations_set%SPECIFICATION(3))
238 SELECT CASE(equations_set_setup%SETUP_TYPE)
244 SELECT CASE(equations_set%SPECIFICATION(3))
249 SELECT CASE(equations_set_setup%ACTION_TYPE)
255 local_error=
"The action type of "//
trim(
number_to_vstring(equations_set_setup%ACTION_TYPE,
"*",err,error))// &
257 &
" is invalid for a standard or quasistatic Darcy equation." 258 CALL flagerror(local_error,err,error,*999)
261 SELECT CASE(equations_set_setup%ACTION_TYPE)
265 equations_set_field_number_of_variables = 1
266 equations_set_field_number_of_components = 2
267 equations_equations_set_field=>equations_set%EQUATIONS_SET_FIELD
268 IF(equations_equations_set_field%EQUATIONS_SET_FIELD_AUTO_CREATED)
THEN 270 CALL field_create_start(equations_set_setup%FIELD_USER_NUMBER,equations_set%REGION, &
271 & equations_equations_set_field%EQUATIONS_SET_FIELD_FIELD,err,error,*999)
272 equations_set_field_field=>equations_equations_set_field%EQUATIONS_SET_FIELD_FIELD
273 CALL field_label_set(equations_set_field_field,
"Equations Set Field",err,error,*999)
274 CALL field_type_set_and_lock(equations_set_field_field,field_general_type,&
276 CALL field_dependent_type_set_and_lock(equations_set_field_field,&
277 & field_independent_type,err,error,*999)
278 CALL field_number_of_variables_set(equations_set_field_field, &
279 & equations_set_field_number_of_variables,err,error,*999)
280 CALL field_variable_types_set_and_lock(equations_set_field_field,&
281 & [field_u_variable_type],err,error,*999)
282 CALL field_dimension_set_and_lock(equations_set_field_field,field_u_variable_type, &
283 & field_vector_dimension_type,err,error,*999)
284 CALL field_data_type_set_and_lock(equations_set_field_field,field_u_variable_type, &
285 & field_intg_type,err,error,*999)
286 CALL field_number_of_components_set_and_lock(equations_set_field_field,&
287 & field_u_variable_type,equations_set_field_number_of_components,err,error,*999)
290 CALL field_type_check(equations_set_setup%FIELD,field_general_type,err,error,*999)
291 CALL field_dependent_type_check(equations_set_setup%FIELD,field_independent_type,err,error,*999)
292 CALL field_number_of_variables_check(equations_set_setup%FIELD,equations_set_field_number_of_variables, &
294 CALL field_variable_types_check(equations_set_setup%FIELD,[field_u_variable_type],err,error,*999)
295 CALL field_dimension_check(equations_set_setup%FIELD,field_u_variable_type,field_vector_dimension_type, &
297 CALL field_data_type_check(equations_set_setup%FIELD,field_u_variable_type,field_intg_type,err,error,*999)
298 CALL field_number_of_components_check(equations_set_setup%FIELD,field_u_variable_type, &
299 & equations_set_field_number_of_components,err,error,*999)
302 IF(equations_set%EQUATIONS_SET_FIELD%EQUATIONS_SET_FIELD_AUTO_CREATED)
THEN 303 CALL field_create_finish(equations_set%EQUATIONS_SET_FIELD%EQUATIONS_SET_FIELD_FIELD,err,error,*999)
304 CALL field_component_values_initialise(equations_set%EQUATIONS_SET_FIELD%EQUATIONS_SET_FIELD_FIELD,&
305 & field_u_variable_type,field_values_set_type, 1, 1_intg, err, error, *999)
306 CALL field_component_values_initialise(equations_set%EQUATIONS_SET_FIELD%EQUATIONS_SET_FIELD_FIELD,&
307 & field_u_variable_type,field_values_set_type, 2, 1_intg, err, error, *999)
311 local_error=
"The action type of "//
trim(
number_to_vstring(equations_set_setup%ACTION_TYPE,
"*",err,error))// &
313 &
" is invalid for a standard or quasistatic Darcy equation." 314 CALL flagerror(local_error,err,error,*999)
322 SELECT CASE(equations_set%SPECIFICATION(3))
331 SELECT CASE(equations_set_setup%ACTION_TYPE)
334 field_variable=>equations_set%GEOMETRY%GEOMETRIC_FIELD%VARIABLE_TYPE_MAP(field_u_variable_type)%PTR
336 CALL field_parametersetensurecreated(equations_set%GEOMETRY%GEOMETRIC_FIELD, field_u_variable_type, &
337 & field_initial_values_set_type, err, error, *999)
339 CALL field_parametersetensurecreated(equations_set%GEOMETRY%GEOMETRIC_FIELD, field_u_variable_type, &
340 & field_previous_values_set_type, err, error, *999)
342 CALL field_parametersetensurecreated(equations_set%GEOMETRY%GEOMETRIC_FIELD, field_u_variable_type, &
343 & field_mesh_displacement_set_type, err, error, *999)
345 CALL field_parametersetensurecreated(equations_set%GEOMETRY%GEOMETRIC_FIELD, field_u_variable_type, &
346 & field_mesh_velocity_set_type, err, error, *999)
348 CALL field_parametersetensurecreated(equations_set%GEOMETRY%GEOMETRIC_FIELD, field_u_variable_type, &
349 & field_negative_mesh_velocity_set_type, err, error, *999)
354 equations_set_field_number_of_components = 2
356 equations_equations_set_field=>equations_set%EQUATIONS_SET_FIELD
357 equations_set_field_field=>equations_equations_set_field%EQUATIONS_SET_FIELD_FIELD
359 IF(equations_equations_set_field%EQUATIONS_SET_FIELD_AUTO_CREATED)
THEN 360 CALL field_mesh_decomposition_get(equations_set%GEOMETRY%GEOMETRIC_FIELD,geometric_decomposition,err,error,*999)
361 CALL field_mesh_decomposition_set_and_lock(equations_set_field_field,&
362 & geometric_decomposition,err,error,*999)
363 CALL field_geometric_field_set_and_lock(equations_set_field_field,&
364 & equations_set%GEOMETRY%GEOMETRIC_FIELD,err,error,*999)
365 CALL field_component_mesh_component_get(equations_set%GEOMETRY%GEOMETRIC_FIELD,field_u_variable_type, &
366 & 1,geometric_component_number,err,error,*999)
367 DO component_idx = 1, equations_set_field_number_of_components
368 CALL field_component_mesh_component_set_and_lock(equations_set%EQUATIONS_SET_FIELD%EQUATIONS_SET_FIELD_FIELD, &
369 & field_u_variable_type,component_idx,geometric_component_number,err,error,*999)
370 CALL field_component_interpolation_set_and_lock(equations_set%EQUATIONS_SET_FIELD%EQUATIONS_SET_FIELD_FIELD, &
371 & field_u_variable_type,component_idx,field_constant_interpolation,err,error,*999)
375 CALL field_scaling_type_get(equations_set%GEOMETRY%GEOMETRIC_FIELD,geometric_scaling_type,err,error,*999)
376 CALL field_scaling_type_set(equations_set%EQUATIONS_SET_FIELD%EQUATIONS_SET_FIELD_FIELD,geometric_scaling_type, &
385 local_error=
"The action type of "//
trim(
number_to_vstring(equations_set_setup%ACTION_TYPE,
"*",err,error))// &
387 &
" is invalid for a linear diffusion equation." 388 CALL flagerror(local_error,err,error,*999)
396 SELECT CASE(equations_set%SPECIFICATION(3))
402 SELECT CASE(equations_set_setup%ACTION_TYPE)
404 IF(equations_set%DEPENDENT%DEPENDENT_FIELD_AUTO_CREATED)
THEN 406 CALL field_create_start(equations_set_setup%FIELD_USER_NUMBER,equations_set%REGION, &
407 & equations_set%DEPENDENT%DEPENDENT_FIELD,err,error,*999)
408 CALL field_type_set_and_lock(equations_set%DEPENDENT%DEPENDENT_FIELD,field_geometric_general_type,err,error,*999)
409 CALL field_dependent_type_set_and_lock(equations_set%DEPENDENT%DEPENDENT_FIELD,field_dependent_type,err,error,*999)
411 CALL field_mesh_decomposition_get(equations_set%GEOMETRY%GEOMETRIC_FIELD,geometric_decomposition,err,error,*999)
412 CALL field_mesh_decomposition_set_and_lock(equations_set%DEPENDENT%DEPENDENT_FIELD,geometric_decomposition, &
414 CALL field_geometric_field_set_and_lock(equations_set%DEPENDENT%DEPENDENT_FIELD, &
415 & equations_set%GEOMETRY%GEOMETRIC_FIELD,err,error,*999)
417 dependent_field_number_of_variables = 2
418 CALL field_number_of_variables_set_and_lock(equations_set%DEPENDENT%DEPENDENT_FIELD, &
419 & dependent_field_number_of_variables,err,error,*999)
420 CALL field_variable_types_set_and_lock(equations_set%DEPENDENT%DEPENDENT_FIELD,[field_u_variable_type, &
421 & field_deludeln_variable_type],err,error,*999)
422 CALL field_variable_label_set(equations_set%DEPENDENT%DEPENDENT_FIELD,field_u_variable_type,
"U",err,error,*999)
423 CALL field_variable_label_set(equations_set%DEPENDENT%DEPENDENT_FIELD,field_deludeln_variable_type,
"del U/del n", &
425 CALL field_dimension_set_and_lock(equations_set%DEPENDENT%DEPENDENT_FIELD,field_u_variable_type, &
426 & field_vector_dimension_type,err,error,*999)
427 CALL field_dimension_set_and_lock(equations_set%DEPENDENT%DEPENDENT_FIELD,field_deludeln_variable_type, &
428 & field_vector_dimension_type,err,error,*999)
429 CALL field_data_type_set_and_lock(equations_set%DEPENDENT%DEPENDENT_FIELD,field_u_variable_type, &
430 & field_dp_type,err,error,*999)
431 CALL field_data_type_set_and_lock(equations_set%DEPENDENT%DEPENDENT_FIELD,field_deludeln_variable_type, &
432 & field_dp_type,err,error,*999)
434 CALL field_number_of_components_get(equations_set%GEOMETRY%GEOMETRIC_FIELD, field_u_variable_type, &
435 & number_of_dimensions, err, error, *999)
436 dependent_field_number_of_components = number_of_dimensions + 1
437 CALL field_number_of_components_set_and_lock(equations_set%DEPENDENT%DEPENDENT_FIELD, field_u_variable_type, &
438 & dependent_field_number_of_components, err, error, *999)
439 CALL field_number_of_components_set_and_lock(equations_set%DEPENDENT%DEPENDENT_FIELD, &
440 & field_deludeln_variable_type, dependent_field_number_of_components, err, error, *999)
442 CALL field_component_mesh_component_get(equations_set%GEOMETRY%GEOMETRIC_FIELD,field_u_variable_type,1, &
443 & geometric_mesh_component,err,error,*999)
444 DO i=1,dependent_field_number_of_components
445 IF( i < dependent_field_number_of_components )
THEN 447 mesh_component = geometric_mesh_component
448 CALL field_component_mesh_component_set(equations_set%DEPENDENT%DEPENDENT_FIELD,field_u_variable_type, &
449 & i, mesh_component,err,error,*999)
450 CALL field_component_mesh_component_set(equations_set%DEPENDENT%DEPENDENT_FIELD,field_deludeln_variable_type, &
451 & i, mesh_component,err,error,*999)
454 mesh_component = geometric_mesh_component
455 CALL field_component_mesh_component_set(equations_set%DEPENDENT%DEPENDENT_FIELD,field_u_variable_type, &
456 & i, mesh_component,err,error,*999)
457 CALL field_component_mesh_component_set(equations_set%DEPENDENT%DEPENDENT_FIELD,field_deludeln_variable_type, &
458 & i, mesh_component,err,error,*999)
462 SELECT CASE(equations_set%SOLUTION_METHOD)
464 DO i = 1, dependent_field_number_of_components
465 CALL field_component_interpolation_set_and_lock(equations_set%DEPENDENT%DEPENDENT_FIELD, &
466 & field_u_variable_type,i,field_node_based_interpolation,err,error,*999)
467 CALL field_component_interpolation_set_and_lock(equations_set%DEPENDENT%DEPENDENT_FIELD, &
468 & field_deludeln_variable_type,i,field_node_based_interpolation,err,error,*999)
471 CALL field_scaling_type_get(equations_set%GEOMETRY%GEOMETRIC_FIELD,geometric_scaling_type,err,error,*999)
472 CALL field_scaling_type_set(equations_set%DEPENDENT%DEPENDENT_FIELD,geometric_scaling_type,err,error,*999)
474 CALL flagerror(
"Not implemented.",err,error,*999)
476 CALL flagerror(
"Not implemented.",err,error,*999)
478 CALL flagerror(
"Not implemented.",err,error,*999)
480 CALL flagerror(
"Not implemented.",err,error,*999)
482 CALL flagerror(
"Not implemented.",err,error,*999)
484 local_error=
"The solution method of "//
trim(
number_to_vstring(equations_set%SOLUTION_METHOD,
"*",err,error))// &
486 CALL flagerror(local_error,err,error,*999)
492 SELECT CASE(equations_set%SPECIFICATION(3))
499 CALL field_type_check(equations_set_setup%FIELD,field_geometric_general_type,err,error,*999)
500 CALL field_dependent_type_check(equations_set_setup%FIELD,field_dependent_type,err,error,*999)
501 CALL field_number_of_variables_check(equations_set_setup%FIELD,4,err,error,*999)
502 CALL field_variable_types_check(equations_set_setup%FIELD,[field_u_variable_type,field_deludeln_variable_type, &
503 & field_v_variable_type,field_delvdeln_variable_type],err,error,*999)
504 CALL field_dimension_check(equations_set_setup%FIELD,field_u_variable_type,field_vector_dimension_type, &
506 CALL field_dimension_check(equations_set_setup%FIELD,field_deludeln_variable_type,field_vector_dimension_type, &
508 CALL field_dimension_check(equations_set_setup%FIELD,field_v_variable_type,field_vector_dimension_type, &
510 CALL field_dimension_check(equations_set_setup%FIELD,field_delvdeln_variable_type,field_vector_dimension_type, &
512 CALL field_data_type_check(equations_set_setup%FIELD,field_u_variable_type,field_dp_type,err,error,*999)
513 CALL field_data_type_check(equations_set_setup%FIELD,field_deludeln_variable_type,field_dp_type,err,error,*999)
514 CALL field_data_type_check(equations_set_setup%FIELD,field_v_variable_type,field_dp_type,err,error,*999)
515 CALL field_data_type_check(equations_set_setup%FIELD,field_delvdeln_variable_type,field_dp_type,err,error,*999)
516 CALL field_number_of_components_get(equations_set%GEOMETRY%GEOMETRIC_FIELD, field_u_variable_type, &
517 & number_of_dimensions, err, error, *999)
519 SELECT CASE(equations_set%SPECIFICATION(3))
521 dependent_field_elasticity_number_of_components = number_of_dimensions
522 dependent_field_darcy_number_of_components = number_of_dimensions + 2
524 dependent_field_elasticity_number_of_components = number_of_dimensions + 1
525 dependent_field_darcy_number_of_components = number_of_dimensions + 1
527 dependent_field_elasticity_number_of_components = number_of_dimensions + 1
528 dependent_field_darcy_number_of_components = number_of_dimensions + 1
531 CALL field_number_of_components_check(equations_set_setup%FIELD,field_u_variable_type, &
532 & dependent_field_elasticity_number_of_components,err,error,*999)
533 CALL field_number_of_components_check(equations_set_setup%FIELD,field_deludeln_variable_type, &
534 & dependent_field_elasticity_number_of_components,err,error,*999)
535 CALL field_number_of_components_check(equations_set_setup%FIELD,field_v_variable_type, &
536 & dependent_field_darcy_number_of_components,err,error,*999)
537 CALL field_number_of_components_check(equations_set_setup%FIELD,field_delvdeln_variable_type, &
538 & dependent_field_darcy_number_of_components,err,error,*999)
540 SELECT CASE(equations_set%SOLUTION_METHOD)
542 CALL field_component_interpolation_check(equations_set_setup%FIELD,field_u_variable_type,1, &
543 & field_node_based_interpolation,err,error,*999)
545 CALL field_component_interpolation_check(equations_set_setup%FIELD,field_deludeln_variable_type,1, &
546 & field_node_based_interpolation,err,error,*999)
547 CALL field_component_interpolation_check(equations_set_setup%FIELD,field_v_variable_type,1, &
548 & field_node_based_interpolation,err,error,*999)
549 CALL field_component_interpolation_check(equations_set_setup%FIELD,field_delvdeln_variable_type,1, &
550 & field_node_based_interpolation,err,error,*999)
552 CALL flagerror(
"Not implemented.",err,error,*999)
554 CALL flagerror(
"Not implemented.",err,error,*999)
556 CALL flagerror(
"Not implemented.",err,error,*999)
558 CALL flagerror(
"Not implemented.",err,error,*999)
560 CALL flagerror(
"Not implemented.",err,error,*999)
562 local_error=
"The solution method of "//
trim(
number_to_vstring(equations_set%SOLUTION_METHOD,
"*",err,error))// &
564 CALL flagerror(local_error,err,error,*999)
571 CALL field_type_check(equations_set_setup%FIELD,field_geometric_general_type,err,error,*999)
572 CALL field_dependent_type_check(equations_set_setup%FIELD,field_dependent_type,err,error,*999)
574 equations_set_field_field=>equations_set%EQUATIONS_SET_FIELD%EQUATIONS_SET_FIELD_FIELD
575 CALL field_parameter_set_data_get(equations_set_field_field,field_u_variable_type, &
576 & field_values_set_type,equations_set_field_data,err,error,*999)
577 ncompartments=equations_set_field_data(2)
578 CALL field_number_of_variables_check(equations_set_setup%FIELD,(2+2*ncompartments),err,error,*999)
579 ALLOCATE(variable_types(2*ncompartments+2))
580 DO num_var=1,ncompartments+1
581 variable_types(2*num_var-1)=field_u_variable_type+(field_number_of_variable_subtypes*(num_var-1))
582 variable_types(2*num_var)=field_deludeln_variable_type+(field_number_of_variable_subtypes*(num_var-1))
584 CALL field_variable_types_check(equations_set_setup%FIELD,variable_types,err,error,*999)
586 CALL field_number_of_components_get(equations_set%GEOMETRY%GEOMETRIC_FIELD,field_u_variable_type, &
587 & number_of_dimensions,err,error,*999)
588 number_of_components=number_of_dimensions+1
589 number_of_darcy_components=number_of_dimensions+1
591 DO num_var=1,2*ncompartments+2
592 CALL field_dimension_check(equations_set_setup%FIELD,variable_types(num_var),field_vector_dimension_type, &
594 CALL field_data_type_check(equations_set_setup%FIELD,variable_types(num_var),field_dp_type,err,error,*999)
595 CALL field_number_of_components_check(equations_set_setup%FIELD,variable_types(num_var),number_of_components, &
599 SELECT CASE(equations_set%SOLUTION_METHOD)
602 DO component_idx=1,number_of_dimensions
603 CALL field_component_interpolation_check(equations_set_setup%FIELD,field_u_variable_type,component_idx, &
604 & field_node_based_interpolation,err,error,*999)
605 CALL field_component_interpolation_check(equations_set_setup%FIELD,field_deludeln_variable_type,component_idx,&
606 & field_node_based_interpolation,err,error,*999)
609 CALL field_component_interpolation_check(equations_set_setup%FIELD,field_u_variable_type,number_of_components,&
610 & field_node_based_interpolation,err,error,*999)
611 CALL field_component_interpolation_check(equations_set_setup%FIELD,field_deludeln_variable_type, &
612 & number_of_components,field_node_based_interpolation,err,error,*999)
613 DO num_var=3,2*ncompartments+2
615 DO component_idx=1,number_of_darcy_components
616 CALL field_component_interpolation_check(equations_set_setup%FIELD,variable_types(num_var),component_idx, &
617 & field_node_based_interpolation,err,error,*999)
621 CALL flagerror(
"Not implemented.",err,error,*999)
623 CALL flagerror(
"Not implemented.",err,error,*999)
625 CALL flagerror(
"Not implemented.",err,error,*999)
627 CALL flagerror(
"Not implemented.",err,error,*999)
629 CALL flagerror(
"Not implemented.",err,error,*999)
631 local_error=
"The solution method of "//
trim(
number_to_vstring(equations_set%SOLUTION_METHOD,
"*",err,error))// &
633 CALL flagerror(local_error,err,error,*999)
637 CALL field_type_check(equations_set_setup%FIELD,field_general_type,err,error,*999)
638 CALL field_dependent_type_check(equations_set_setup%FIELD,field_dependent_type,err,error,*999)
639 equations_set_field_field=>equations_set%EQUATIONS_SET_FIELD%EQUATIONS_SET_FIELD_FIELD
640 CALL field_parameter_set_data_get(equations_set_field_field,field_u_variable_type, &
641 & field_values_set_type,equations_set_field_data,err,error,*999)
642 ncompartments=equations_set_field_data(2)
643 CALL field_number_of_variables_check(equations_set_setup%FIELD,2*ncompartments,err,error,*999)
645 ALLOCATE(variable_types(2*ncompartments))
646 DO num_var=1,ncompartments
647 variable_types(2*num_var-1)=field_u_variable_type+(field_number_of_variable_subtypes*(num_var-1))
648 variable_types(2*num_var)=field_deludeln_variable_type+(field_number_of_variable_subtypes*(num_var-1))
650 CALL field_variable_types_check(equations_set_setup%FIELD,variable_types,err,error,*999)
652 CALL field_number_of_components_get(equations_set%GEOMETRY%GEOMETRIC_FIELD,field_u_variable_type, &
653 & number_of_dimensions,err,error,*999)
654 DO num_var=1,2*ncompartments
655 CALL field_dimension_check(equations_set_setup%FIELD,variable_types(num_var), &
656 & field_vector_dimension_type,err,error,*999)
657 CALL field_data_type_check(equations_set_setup%FIELD,variable_types(num_var),field_dp_type,err,error,*999)
658 CALL field_number_of_components_check(equations_set_setup%FIELD,variable_types(num_var), &
659 & number_of_dimensions+1,err,error,*999)
662 SELECT CASE(equations_set%SOLUTION_METHOD)
665 DO num_var=1,2*ncompartments
666 DO component_idx=1,number_of_dimensions+1
667 CALL field_component_interpolation_check(equations_set_setup%FIELD,variable_types(num_var),component_idx, &
668 & field_node_based_interpolation,err,error,*999)
673 CALL flagerror(
"Not implemented.",err,error,*999)
675 CALL flagerror(
"Not implemented.",err,error,*999)
677 CALL flagerror(
"Not implemented.",err,error,*999)
679 CALL flagerror(
"Not implemented.",err,error,*999)
681 CALL flagerror(
"Not implemented.",err,error,*999)
683 local_error=
"The solution method of "//
trim(
number_to_vstring(equations_set%SOLUTION_METHOD,
"*",err,error))// &
685 CALL flagerror(local_error,err,error,*999)
692 CALL field_type_check(equations_set_setup%FIELD,field_general_type,err,error,*999)
693 CALL field_dependent_type_check(equations_set_setup%FIELD,field_dependent_type,err,error,*999)
694 CALL field_number_of_variables_check(equations_set_setup%FIELD,2,err,error,*999)
695 CALL field_variable_types_check(equations_set_setup%FIELD,[field_u_variable_type,field_deludeln_variable_type],&
697 CALL field_dimension_check(equations_set_setup%FIELD,field_u_variable_type,field_vector_dimension_type, &
699 CALL field_dimension_check(equations_set_setup%FIELD,field_deludeln_variable_type,field_vector_dimension_type, &
701 CALL field_data_type_check(equations_set_setup%FIELD,field_u_variable_type,field_dp_type,err,error,*999)
702 CALL field_data_type_check(equations_set_setup%FIELD,field_deludeln_variable_type,field_dp_type,err,error,*999)
703 CALL field_number_of_components_get(equations_set%GEOMETRY%GEOMETRIC_FIELD, field_u_variable_type, &
704 & number_of_dimensions, err, error, *999)
705 dependent_field_number_of_components = number_of_dimensions + 1
706 CALL field_number_of_components_check(equations_set_setup%FIELD,field_u_variable_type, &
707 & dependent_field_number_of_components,err,error,*999)
708 CALL field_number_of_components_check(equations_set_setup%FIELD,field_deludeln_variable_type, &
709 & dependent_field_number_of_components,err,error,*999)
711 SELECT CASE(equations_set%SOLUTION_METHOD)
713 CALL field_component_interpolation_check(equations_set_setup%FIELD,field_u_variable_type,1, &
714 & field_node_based_interpolation,err,error,*999)
715 CALL field_component_interpolation_check(equations_set_setup%FIELD,field_deludeln_variable_type,1, &
716 & field_node_based_interpolation,err,error,*999)
718 CALL flagerror(
"Not implemented.",err,error,*999)
720 CALL flagerror(
"Not implemented.",err,error,*999)
722 CALL flagerror(
"Not implemented.",err,error,*999)
724 CALL flagerror(
"Not implemented.",err,error,*999)
726 CALL flagerror(
"Not implemented.",err,error,*999)
728 local_error=
"The solution method of "//
trim(
number_to_vstring(equations_set%SOLUTION_METHOD,
"*",err,error))// &
730 CALL flagerror(local_error,err,error,*999)
735 IF(equations_set%DEPENDENT%DEPENDENT_FIELD_AUTO_CREATED)
THEN 736 CALL field_create_finish(equations_set%DEPENDENT%DEPENDENT_FIELD,err,error,*999)
737 CALL field_parameter_set_create(equations_set%DEPENDENT%DEPENDENT_FIELD,field_u_variable_type, &
738 & field_initial_values_set_type,err,error,*999)
742 CALL field_parameter_set_create(equations_set%DEPENDENT%DEPENDENT_FIELD,field_u_variable_type, &
743 & field_relative_velocity_set_type,err,error,*999)
746 local_error=
"The action type of "//
trim(
number_to_vstring(equations_set_setup%ACTION_TYPE,
"*",err,error))// &
748 &
" is invalid for a standard, quasistatic or ALE Darcy equation" 749 CALL flagerror(local_error,err,error,*999)
757 SELECT CASE(equations_set%SPECIFICATION(3))
764 SELECT CASE(equations_set_setup%ACTION_TYPE)
766 IF(equations_set%INDEPENDENT%INDEPENDENT_FIELD_AUTO_CREATED)
THEN 768 CALL field_create_start(equations_set_setup%FIELD_USER_NUMBER,equations_set%REGION, &
769 & equations_set%INDEPENDENT%INDEPENDENT_FIELD,err,error,*999)
770 CALL field_type_set_and_lock(equations_set%INDEPENDENT%INDEPENDENT_FIELD,field_general_type,err,error,*999)
772 CALL field_dependent_type_set_and_lock(equations_set%INDEPENDENT%INDEPENDENT_FIELD,field_independent_type,&
775 CALL field_mesh_decomposition_get(equations_set%GEOMETRY%GEOMETRIC_FIELD,geometric_decomposition,err,error,*999)
776 CALL field_mesh_decomposition_set_and_lock(equations_set%INDEPENDENT%INDEPENDENT_FIELD,geometric_decomposition, &
778 CALL field_geometric_field_set_and_lock(equations_set%INDEPENDENT%INDEPENDENT_FIELD, &
779 & equations_set%GEOMETRY%GEOMETRIC_FIELD,err,error,*999)
781 independent_field_number_of_variables = 2
782 CALL field_number_of_variables_set_and_lock(equations_set%INDEPENDENT%INDEPENDENT_FIELD, &
783 & independent_field_number_of_variables,err,error,*999)
784 CALL field_variable_types_set_and_lock(equations_set%INDEPENDENT%INDEPENDENT_FIELD,[field_u_variable_type, &
785 & field_deludeln_variable_type],err,error,*999)
786 CALL field_variable_label_set(equations_set%INDEPENDENT%INDEPENDENT_FIELD,field_u_variable_type,
"Independent U", &
788 CALL field_variable_label_set(equations_set%INDEPENDENT%INDEPENDENT_FIELD,field_deludeln_variable_type, &
789 &
"Independent del U/del n",err,error,*999)
790 CALL field_dimension_set_and_lock(equations_set%INDEPENDENT%INDEPENDENT_FIELD,field_u_variable_type, &
791 & field_vector_dimension_type,err,error,*999)
792 CALL field_dimension_set_and_lock(equations_set%INDEPENDENT%INDEPENDENT_FIELD,field_deludeln_variable_type, &
793 & field_vector_dimension_type,err,error,*999)
794 CALL field_data_type_set_and_lock(equations_set%INDEPENDENT%INDEPENDENT_FIELD,field_u_variable_type, &
795 & field_dp_type,err,error,*999)
796 CALL field_data_type_set_and_lock(equations_set%INDEPENDENT%INDEPENDENT_FIELD,field_deludeln_variable_type, &
797 & field_dp_type,err,error,*999)
799 CALL field_number_of_components_get(equations_set%GEOMETRY%GEOMETRIC_FIELD, field_u_variable_type, &
800 & number_of_dimensions, err, error, *999)
801 independent_field_number_of_components = number_of_dimensions
802 CALL field_number_of_components_set_and_lock(equations_set%INDEPENDENT%INDEPENDENT_FIELD, field_u_variable_type, &
803 & independent_field_number_of_components, err, error, *999)
804 CALL field_number_of_components_set_and_lock(equations_set%INDEPENDENT%INDEPENDENT_FIELD, &
805 & field_deludeln_variable_type, independent_field_number_of_components, err, error, *999)
807 CALL field_component_mesh_component_get(equations_set%GEOMETRY%GEOMETRIC_FIELD,field_u_variable_type,1, &
808 & geometric_mesh_component,err,error,*999)
809 DO i=1,independent_field_number_of_components
810 IF( i < independent_field_number_of_components )
THEN 812 mesh_component = geometric_mesh_component
813 CALL field_component_mesh_component_set(equations_set%INDEPENDENT%INDEPENDENT_FIELD,field_u_variable_type, &
814 & i, mesh_component,err,error,*999)
815 CALL field_component_mesh_component_set(equations_set%INDEPENDENT%INDEPENDENT_FIELD,&
816 & field_deludeln_variable_type, i, mesh_component,err,error,*999)
819 mesh_component = geometric_mesh_component
820 CALL field_component_mesh_component_set(equations_set%INDEPENDENT%INDEPENDENT_FIELD,field_u_variable_type, &
821 & i, mesh_component,err,error,*999)
822 CALL field_component_mesh_component_set(equations_set%INDEPENDENT%INDEPENDENT_FIELD,&
823 & field_deludeln_variable_type, i, mesh_component,err,error,*999)
827 SELECT CASE(equations_set%SOLUTION_METHOD)
829 DO i = 1, independent_field_number_of_components
830 CALL field_component_interpolation_set_and_lock(equations_set%INDEPENDENT%INDEPENDENT_FIELD, &
831 & field_u_variable_type,i,field_node_based_interpolation,err,error,*999)
832 CALL field_component_interpolation_set_and_lock(equations_set%INDEPENDENT%INDEPENDENT_FIELD, &
833 & field_deludeln_variable_type,i,field_node_based_interpolation,err,error,*999)
836 CALL field_scaling_type_get(equations_set%GEOMETRY%GEOMETRIC_FIELD,geometric_scaling_type,err,error,*999)
837 CALL field_scaling_type_set(equations_set%INDEPENDENT%INDEPENDENT_FIELD,geometric_scaling_type,err,error,*999)
839 CALL flagerror(
"Not implemented.",err,error,*999)
841 CALL flagerror(
"Not implemented.",err,error,*999)
843 CALL flagerror(
"Not implemented.",err,error,*999)
845 CALL flagerror(
"Not implemented.",err,error,*999)
847 CALL flagerror(
"Not implemented.",err,error,*999)
849 local_error=
"The solution method of "//
trim(
number_to_vstring(equations_set%SOLUTION_METHOD,
"*",err,error))// &
851 CALL flagerror(local_error,err,error,*999)
855 CALL field_type_check(equations_set_setup%FIELD,field_general_type,err,error,*999)
857 CALL field_dependent_type_check(equations_set_setup%FIELD,field_independent_type,err,error,*999)
859 CALL field_number_of_variables_check(equations_set_setup%FIELD,2,err,error,*999)
860 CALL field_variable_types_check(equations_set_setup%FIELD,[field_u_variable_type,field_deludeln_variable_type], &
862 CALL field_dimension_check(equations_set_setup%FIELD,field_u_variable_type,field_vector_dimension_type, &
864 CALL field_dimension_check(equations_set_setup%FIELD,field_deludeln_variable_type,field_vector_dimension_type, &
866 CALL field_data_type_check(equations_set_setup%FIELD,field_u_variable_type,field_dp_type,err,error,*999)
867 CALL field_data_type_check(equations_set_setup%FIELD,field_deludeln_variable_type,field_dp_type,err,error,*999)
868 CALL field_number_of_components_get(equations_set%GEOMETRY%GEOMETRIC_FIELD, field_u_variable_type, &
869 & number_of_dimensions, err, error, *999)
870 independent_field_number_of_components = number_of_dimensions
871 CALL field_number_of_components_check(equations_set_setup%FIELD,field_u_variable_type, &
872 & independent_field_number_of_components,err,error,*999)
873 CALL field_number_of_components_check(equations_set_setup%FIELD,field_deludeln_variable_type, &
874 & independent_field_number_of_components,err,error,*999)
876 SELECT CASE(equations_set%SOLUTION_METHOD)
878 CALL field_component_interpolation_check(equations_set_setup%FIELD,field_u_variable_type,1, &
879 & field_node_based_interpolation,err,error,*999)
880 CALL field_component_interpolation_check(equations_set_setup%FIELD,field_deludeln_variable_type,1, &
881 & field_node_based_interpolation,err,error,*999)
883 CALL flagerror(
"Not implemented.",err,error,*999)
885 CALL flagerror(
"Not implemented.",err,error,*999)
887 CALL flagerror(
"Not implemented.",err,error,*999)
889 CALL flagerror(
"Not implemented.",err,error,*999)
891 CALL flagerror(
"Not implemented.",err,error,*999)
893 local_error=
"The solution method of "//
trim(
number_to_vstring(equations_set%SOLUTION_METHOD,
"*",err,error))// &
895 CALL flagerror(local_error,err,error,*999)
899 IF(equations_set%INDEPENDENT%INDEPENDENT_FIELD_AUTO_CREATED)
THEN 900 CALL field_create_finish(equations_set%INDEPENDENT%INDEPENDENT_FIELD,err,error,*999)
903 local_error=
"The action type of "//
trim(
number_to_vstring(equations_set_setup%ACTION_TYPE,
"*",err,error))// &
905 &
" is invalid for a standard, quasistatic or ALE Darcy equation" 906 CALL flagerror(local_error,err,error,*999)
914 SELECT CASE(equations_set%SPECIFICATION(3))
919 material_field_number_of_variables = 1
920 SELECT CASE(equations_set%SPECIFICATION(3))
923 material_field_number_of_components = 2
926 material_field_number_of_components = 7
928 SELECT CASE(equations_set_setup%ACTION_TYPE)
930 equations_materials=>equations_set%MATERIALS
931 IF(
ASSOCIATED(equations_materials))
THEN 932 IF(equations_materials%MATERIALS_FIELD_AUTO_CREATED)
THEN 934 CALL field_create_start(equations_set_setup%FIELD_USER_NUMBER,equations_set%REGION,equations_materials% &
935 & materials_field,err,error,*999)
936 CALL field_type_set_and_lock(equations_materials%MATERIALS_FIELD,field_material_type,err,error,*999)
937 CALL field_dependent_type_set_and_lock(equations_materials%MATERIALS_FIELD,field_independent_type,err,error,*999)
938 CALL field_mesh_decomposition_get(equations_set%GEOMETRY%GEOMETRIC_FIELD,geometric_decomposition,err,error,*999)
939 CALL field_mesh_decomposition_set_and_lock(equations_materials%MATERIALS_FIELD,geometric_decomposition, &
941 CALL field_geometric_field_set_and_lock(equations_materials%MATERIALS_FIELD,equations_set%GEOMETRY% &
942 & geometric_field,err,error,*999)
943 CALL field_number_of_variables_set(equations_materials%MATERIALS_FIELD, &
944 & material_field_number_of_variables,err,error,*999)
945 CALL field_variable_types_set_and_lock(equations_materials%MATERIALS_FIELD,[field_u_variable_type], &
947 CALL field_variable_label_set(equations_materials%MATERIALS_FIELD,field_u_variable_type,
"Material", &
949 CALL field_dimension_set_and_lock(equations_materials%MATERIALS_FIELD,field_u_variable_type, &
950 & field_vector_dimension_type,err,error,*999)
951 CALL field_data_type_set_and_lock(equations_materials%MATERIALS_FIELD,field_u_variable_type, &
952 & field_dp_type,err,error,*999)
953 CALL field_number_of_components_set_and_lock(equations_materials%MATERIALS_FIELD,field_u_variable_type, &
954 & material_field_number_of_components,err,error,*999)
955 CALL field_component_mesh_component_get(equations_set%GEOMETRY%GEOMETRIC_FIELD,field_u_variable_type, &
956 & 1,geometric_component_number,err,error,*999)
960 DO i = 1, material_field_number_of_components
961 CALL field_component_interpolation_set(equations_materials%MATERIALS_FIELD,field_u_variable_type, &
962 & i,field_node_based_interpolation,err,error,*999)
963 CALL field_component_mesh_component_set(equations_materials%MATERIALS_FIELD,field_u_variable_type, &
964 & i,geometric_component_number,err,error,*999)
968 CALL field_scaling_type_get(equations_set%GEOMETRY%GEOMETRIC_FIELD,geometric_scaling_type,err,error,*999)
969 CALL field_scaling_type_set(equations_materials%MATERIALS_FIELD,geometric_scaling_type,err,error,*999)
972 CALL field_type_check(equations_set_setup%FIELD,field_material_type,err,error,*999)
973 CALL field_dependent_type_check(equations_set_setup%FIELD,field_independent_type,err,error,*999)
974 CALL field_number_of_variables_check(equations_set_setup%FIELD,material_field_number_of_variables,err,error,*999)
975 CALL field_variable_types_check(equations_set_setup%FIELD,[field_u_variable_type],err,error,*999)
976 CALL field_dimension_check(equations_set_setup%FIELD,field_u_variable_type,field_vector_dimension_type, &
978 CALL field_data_type_check(equations_set_setup%FIELD,field_u_variable_type,field_dp_type,err,error,*999)
979 CALL field_number_of_components_check(equations_set_setup%FIELD,field_u_variable_type, &
980 & material_field_number_of_components,err,error,*999)
983 CALL flagerror(
"Equations set materials is not associated.",err,error,*999)
986 equations_materials=>equations_set%MATERIALS
987 IF(
ASSOCIATED(equations_materials) )
THEN 988 IF( equations_materials%MATERIALS_FIELD_AUTO_CREATED )
THEN 989 CALL field_create_finish(equations_materials%MATERIALS_FIELD,err,error,*999)
991 DO i=1,material_field_number_of_components
992 CALL field_component_values_initialise(equations_materials%MATERIALS_FIELD,field_u_variable_type, &
993 & field_values_set_type, i, 1.0_dp, err, error, *999)
997 CALL flagerror(
"Equations set materials is not associated.",err,error,*999)
1000 local_error=
"The action type of "//
trim(
number_to_vstring(equations_set_setup%ACTION_TYPE,
"*",err,error))// &
1001 &
" for a setup type of "//
trim(
number_to_vstring(equations_set_setup%SETUP_TYPE,
"*",err,error))// &
1002 &
" is invalid for a standard, quasistatic or ALE Darcy equation." 1003 CALL flagerror(local_error,err,error,*999)
1010 equations_set_field_field=>equations_set%EQUATIONS_SET_FIELD%EQUATIONS_SET_FIELD_FIELD
1011 CALL field_parameter_set_data_get(equations_set_field_field,field_u_variable_type, &
1012 & field_values_set_type,equations_set_field_data,err,error,*999)
1013 ncompartments=equations_set_field_data(2)
1014 material_field_number_of_variables = 3
1015 material_field_number_of_u_var_components = 2
1016 material_field_number_of_v_var_components = ncompartments
1017 material_field_number_of_u1_var_components = 3
1018 SELECT CASE(equations_set_setup%ACTION_TYPE)
1020 equations_materials=>equations_set%MATERIALS
1021 IF(
ASSOCIATED(equations_materials))
THEN 1022 IF(equations_materials%MATERIALS_FIELD_AUTO_CREATED)
THEN 1024 CALL field_create_start(equations_set_setup%FIELD_USER_NUMBER,equations_set%REGION,equations_materials% &
1025 & materials_field,err,error,*999)
1026 CALL field_type_set_and_lock(equations_materials%MATERIALS_FIELD,field_material_type,err,error,*999)
1027 CALL field_dependent_type_set_and_lock(equations_materials%MATERIALS_FIELD,field_independent_type,err,error,*999)
1028 CALL field_mesh_decomposition_get(equations_set%GEOMETRY%GEOMETRIC_FIELD,geometric_decomposition,err,error,*999)
1029 CALL field_mesh_decomposition_set_and_lock(equations_materials%MATERIALS_FIELD,geometric_decomposition, &
1031 CALL field_geometric_field_set_and_lock(equations_materials%MATERIALS_FIELD,equations_set%GEOMETRY% &
1032 & geometric_field,err,error,*999)
1033 CALL field_number_of_variables_set(equations_materials%MATERIALS_FIELD, &
1034 & material_field_number_of_variables,err,error,*999)
1035 CALL field_variable_types_set_and_lock(equations_materials%MATERIALS_FIELD,[field_u_variable_type, &
1036 & field_v_variable_type,field_u1_variable_type],err,error,*999)
1037 CALL field_dimension_set_and_lock(equations_materials%MATERIALS_FIELD,field_u_variable_type, &
1038 & field_vector_dimension_type,err,error,*999)
1039 CALL field_data_type_set_and_lock(equations_materials%MATERIALS_FIELD,field_u_variable_type, &
1040 & field_dp_type,err,error,*999)
1041 CALL field_dimension_set_and_lock(equations_materials%MATERIALS_FIELD,field_v_variable_type, &
1042 & field_vector_dimension_type,err,error,*999)
1043 CALL field_data_type_set_and_lock(equations_materials%MATERIALS_FIELD,field_v_variable_type, &
1044 & field_dp_type,err,error,*999)
1045 CALL field_dimension_set_and_lock(equations_materials%MATERIALS_FIELD,field_u1_variable_type, &
1046 & field_vector_dimension_type,err,error,*999)
1047 CALL field_data_type_set_and_lock(equations_materials%MATERIALS_FIELD,field_u1_variable_type, &
1048 & field_dp_type,err,error,*999)
1049 CALL field_number_of_components_set_and_lock(equations_materials%MATERIALS_FIELD,field_u_variable_type, &
1050 & material_field_number_of_u_var_components,err,error,*999)
1051 CALL field_number_of_components_set_and_lock(equations_materials%MATERIALS_FIELD,field_v_variable_type, &
1052 & material_field_number_of_v_var_components,err,error,*999)
1053 CALL field_number_of_components_set_and_lock(equations_materials%MATERIALS_FIELD,field_u1_variable_type, &
1054 & material_field_number_of_u1_var_components,err,error,*999)
1055 CALL field_component_mesh_component_get(equations_set%GEOMETRY%GEOMETRIC_FIELD,field_u_variable_type, &
1056 & 1,geometric_component_number,err,error,*999)
1060 DO i = 1, material_field_number_of_u_var_components
1061 CALL field_component_interpolation_set(equations_materials%MATERIALS_FIELD,field_u_variable_type, &
1062 & i,field_node_based_interpolation,err,error,*999)
1063 CALL field_component_mesh_component_set(equations_materials%MATERIALS_FIELD,field_u_variable_type, &
1064 & i,geometric_component_number,err,error,*999)
1066 DO i = 1, material_field_number_of_v_var_components
1067 CALL field_component_interpolation_set(equations_materials%MATERIALS_FIELD,field_v_variable_type, &
1068 & i,field_node_based_interpolation,err,error,*999)
1069 CALL field_component_mesh_component_set(equations_materials%MATERIALS_FIELD,field_v_variable_type, &
1070 & i,geometric_component_number,err,error,*999)
1072 DO i = 1, material_field_number_of_u1_var_components
1073 CALL field_component_interpolation_set(equations_materials%MATERIALS_FIELD,field_u1_variable_type, &
1074 & i,field_node_based_interpolation,err,error,*999)
1075 CALL field_component_mesh_component_set(equations_materials%MATERIALS_FIELD,field_u1_variable_type, &
1076 & i,geometric_component_number,err,error,*999)
1080 CALL field_scaling_type_get(equations_set%GEOMETRY%GEOMETRIC_FIELD,geometric_scaling_type,err,error,*999)
1081 CALL field_scaling_type_set(equations_materials%MATERIALS_FIELD,geometric_scaling_type,err,error,*999)
1084 CALL field_type_check(equations_set_setup%FIELD,field_material_type,err,error,*999)
1085 CALL field_dependent_type_check(equations_set_setup%FIELD,field_independent_type,err,error,*999)
1086 CALL field_number_of_variables_check(equations_set_setup%FIELD,material_field_number_of_variables,err,error,*999)
1087 CALL field_variable_types_check(equations_set_setup%FIELD,[field_u_variable_type, &
1088 & field_v_variable_type,field_u1_variable_type],err,error,*999)
1089 CALL field_dimension_check(equations_set_setup%FIELD,field_u_variable_type,field_vector_dimension_type, &
1091 CALL field_dimension_check(equations_set_setup%FIELD,field_v_variable_type,field_vector_dimension_type, &
1093 CALL field_dimension_check(equations_set_setup%FIELD,field_u1_variable_type,field_vector_dimension_type, &
1095 CALL field_data_type_check(equations_set_setup%FIELD,field_u_variable_type,field_dp_type,err,error,*999)
1096 CALL field_data_type_check(equations_set_setup%FIELD,field_v_variable_type,field_dp_type,err,error,*999)
1097 CALL field_data_type_check(equations_set_setup%FIELD,field_u1_variable_type,field_dp_type,err,error,*999)
1098 CALL field_number_of_components_check(equations_set_setup%FIELD,field_u_variable_type, &
1099 & material_field_number_of_u_var_components,err,error,*999)
1100 CALL field_number_of_components_check(equations_set_setup%FIELD,field_v_variable_type, &
1101 & material_field_number_of_v_var_components,err,error,*999)
1102 CALL field_number_of_components_check(equations_set_setup%FIELD,field_u1_variable_type, &
1103 & material_field_number_of_u1_var_components,err,error,*999)
1106 CALL flagerror(
"Equations set materials is not associated.",err,error,*999)
1109 equations_materials=>equations_set%MATERIALS
1110 IF(
ASSOCIATED(equations_materials) )
THEN 1111 IF( equations_materials%MATERIALS_FIELD_AUTO_CREATED )
THEN 1112 CALL field_create_finish(equations_materials%MATERIALS_FIELD,err,error,*999)
1114 DO i=1,material_field_number_of_u_var_components
1115 CALL field_component_values_initialise(equations_materials%MATERIALS_FIELD,field_u_variable_type, &
1116 & field_values_set_type, i, 1.0_dp, err, error, *999)
1118 DO i=1,material_field_number_of_v_var_components
1119 CALL field_component_values_initialise(equations_materials%MATERIALS_FIELD,field_v_variable_type, &
1120 & field_values_set_type, i, 0.0_dp, err, error, *999)
1122 DO i=1,material_field_number_of_u1_var_components
1123 CALL field_component_values_initialise(equations_materials%MATERIALS_FIELD,field_u1_variable_type, &
1124 & field_values_set_type, i, 0.0_dp, err, error, *999)
1128 CALL flagerror(
"Equations set materials is not associated.",err,error,*999)
1131 local_error=
"The action type of "//
trim(
number_to_vstring(equations_set_setup%ACTION_TYPE,
"*",err,error))// &
1132 &
" for a setup type of "//
trim(
number_to_vstring(equations_set_setup%SETUP_TYPE,
"*",err,error))// &
1133 &
" is invalid for a standard, quasistatic or ALE Darcy equation." 1134 CALL flagerror(local_error,err,error,*999)
1144 SELECT CASE(equations_set%SPECIFICATION(3))
1146 SELECT CASE(equations_set_setup%ACTION_TYPE)
1149 IF(equations_set%DEPENDENT%DEPENDENT_FINISHED)
THEN 1150 IF(
ASSOCIATED(equations_set%DEPENDENT%DEPENDENT_FIELD))
THEN 1151 IF(
ASSOCIATED(equations_set%GEOMETRY%GEOMETRIC_FIELD))
THEN 1152 CALL field_number_of_components_get(equations_set%GEOMETRY%GEOMETRIC_FIELD,field_u_variable_type, &
1153 & number_of_dimensions,err,error,*999)
1154 SELECT CASE(equations_set_setup%ANALYTIC_FUNCTION_TYPE)
1177 local_error=
"The specified analytic function type of "// &
1179 &
" is invalid for an analytic Darcy problem." 1180 CALL flagerror(local_error,err,error,*999)
1183 CALL flagerror(
"Equations set geometric field is not associated.",err,error,*999)
1186 CALL flagerror(
"Equations set dependent field is not associated.",err,error,*999)
1189 CALL flagerror(
"Equations set dependent field has not been finished.",err,error,*999)
1192 IF(
ASSOCIATED(equations_set%ANALYTIC))
THEN 1193 IF(
ASSOCIATED(equations_set%ANALYTIC%ANALYTIC_FIELD))
THEN 1194 IF(equations_set%ANALYTIC%ANALYTIC_FIELD_AUTO_CREATED)
THEN 1196 CALL field_create_finish(equations_set%DEPENDENT%DEPENDENT_FIELD,err,error,*999)
1200 CALL flagerror(
"Equations set analytic is not associated.",err,error,*999)
1203 local_error=
"The action type of "//
trim(
number_to_vstring(equations_set_setup%ACTION_TYPE,
"*",err,error))// &
1204 &
" for a setup type of "//
trim(
number_to_vstring(equations_set_setup%SETUP_TYPE,
"*",err,error))// &
1205 &
" is invalid for an analytic Darcy problem." 1206 CALL flagerror(local_error,err,error,*999)
1209 SELECT CASE(equations_set_setup%ACTION_TYPE)
1212 IF(equations_set%DEPENDENT%DEPENDENT_FINISHED)
THEN 1213 IF(
ASSOCIATED(equations_set%DEPENDENT%DEPENDENT_FIELD))
THEN 1214 IF(
ASSOCIATED(equations_set%GEOMETRY%GEOMETRIC_FIELD))
THEN 1215 CALL field_number_of_components_get(equations_set%GEOMETRY%GEOMETRIC_FIELD,field_u_variable_type, &
1216 & number_of_dimensions,err,error,*999)
1218 equations_set%ANALYTIC%ANALYTIC_USER_PARAMS(1)=0.0_dp
1219 SELECT CASE(equations_set_setup%ANALYTIC_FUNCTION_TYPE)
1224 local_error=
"The specified analytic function type of "// &
1226 &
" is invalid for an analytic Darcy problem." 1227 CALL flagerror(local_error,err,error,*999)
1230 CALL flagerror(
"Equations set geometric field is not associated.",err,error,*999)
1233 CALL flagerror(
"Equations set dependent field is not associated.",err,error,*999)
1236 CALL flagerror(
"Equations set dependent field has not been finished.",err,error,*999)
1239 IF(
ASSOCIATED(equations_set%ANALYTIC))
THEN 1240 IF(
ASSOCIATED(equations_set%ANALYTIC%ANALYTIC_FIELD))
THEN 1241 IF(equations_set%ANALYTIC%ANALYTIC_FIELD_AUTO_CREATED)
THEN 1243 CALL field_create_finish(equations_set%DEPENDENT%DEPENDENT_FIELD,err,error,*999)
1247 CALL flagerror(
"Equations set analytic is not associated.",err,error,*999)
1250 local_error=
"The action type of "//
trim(
number_to_vstring(equations_set_setup%ACTION_TYPE,
"*",err,error))// &
1251 &
" for a setup type of "//
trim(
number_to_vstring(equations_set_setup%SETUP_TYPE,
"*",err,error))// &
1252 &
" is invalid for an analytic Darcy problem." 1253 CALL flagerror(local_error,err,error,*999)
1256 local_error=
"The equation set subtype of "//
trim(
number_to_vstring(equations_set%SPECIFICATION(3),
"*",err,error))// &
1257 &
" for a setup type of "//
trim(
number_to_vstring(equations_set_setup%SETUP_TYPE,
"*",err,error))// &
1258 &
" is invalid for a Darcy equation." 1259 CALL flagerror(local_error,err,error,*999)
1266 SELECT CASE(equations_set%SPECIFICATION(3))
1272 SELECT CASE(equations_set_setup%ACTION_TYPE)
1274 equations_source=>equations_set%SOURCE
1275 IF(
ASSOCIATED(equations_source))
THEN 1276 IF(equations_source%SOURCE_FIELD_AUTO_CREATED)
THEN 1277 CALL field_create_start(equations_set_setup%FIELD_USER_NUMBER,equations_set%REGION,equations_source% &
1279 CALL field_label_set(equations_source%SOURCE_FIELD,
"Source Field",err,error,*999)
1280 CALL field_type_set_and_lock(equations_source%SOURCE_FIELD,field_general_type,err,error,*999)
1281 CALL field_dependent_type_set_and_lock(equations_source%SOURCE_FIELD,field_independent_type,err,error,*999)
1282 CALL field_mesh_decomposition_get(equations_set%GEOMETRY%GEOMETRIC_FIELD,geometric_decomposition,err,error,*999)
1283 CALL field_mesh_decomposition_set_and_lock(equations_source%SOURCE_FIELD,geometric_decomposition, &
1285 CALL field_geometric_field_set_and_lock(equations_source%SOURCE_FIELD,equations_set%GEOMETRY% &
1286 & geometric_field,err,error,*999)
1287 CALL field_number_of_variables_set_and_lock(equations_source%SOURCE_FIELD,1,err,error,*999)
1288 CALL field_variable_types_set_and_lock(equations_source%SOURCE_FIELD,[field_u_variable_type], &
1290 CALL field_variable_label_set(equations_source%SOURCE_FIELD,field_u_variable_type,
"Source", &
1292 CALL field_dimension_set_and_lock(equations_source%SOURCE_FIELD,field_u_variable_type, &
1293 & field_vector_dimension_type,err,error,*999)
1294 CALL field_data_type_set_and_lock(equations_source%SOURCE_FIELD,field_u_variable_type, &
1295 & field_dp_type,err,error,*999)
1296 CALL field_number_of_components_get(equations_set%GEOMETRY%GEOMETRIC_FIELD, field_u_variable_type, &
1297 & number_of_dimensions, err, error, *999)
1298 number_of_source_components = number_of_dimensions + 1
1299 CALL field_number_of_components_set_and_lock(equations_source%SOURCE_FIELD,field_u_variable_type, &
1300 & number_of_source_components,err,error,*999)
1306 DO component_idx=1,number_of_dimensions
1307 CALL field_component_mesh_component_get(equations_set%GEOMETRY%GEOMETRIC_FIELD,field_u_variable_type, &
1308 & component_idx,geometric_mesh_component,err,error,*999)
1309 CALL field_component_interpolation_set(equations_source%SOURCE_FIELD,field_u_variable_type, &
1310 & component_idx,field_node_based_interpolation,err,error,*999)
1311 CALL field_component_mesh_component_set(equations_source%SOURCE_FIELD,field_u_variable_type, &
1312 & component_idx,geometric_mesh_component,err,error,*999)
1315 CALL field_component_interpolation_set(equations_source%SOURCE_FIELD,field_u_variable_type, &
1316 & number_of_dimensions + 1,field_node_based_interpolation,err,error,*999)
1317 CALL field_component_mesh_component_set(equations_source%SOURCE_FIELD,field_u_variable_type, &
1318 & number_of_dimensions + 1,geometric_mesh_component,err,error,*999)
1321 CALL field_scaling_type_get(equations_set%GEOMETRY%GEOMETRIC_FIELD,geometric_scaling_type,err,error,*999)
1322 CALL field_scaling_type_set(equations_source%SOURCE_FIELD,geometric_scaling_type,err,error,*999)
1325 CALL field_type_check(equations_set_setup%FIELD,field_general_type,err,error,*999)
1326 CALL field_dependent_type_check(equations_set_setup%FIELD,field_independent_type,err,error,*999)
1327 CALL field_number_of_variables_check(equations_set_setup%FIELD,1,err,error,*999)
1328 CALL field_variable_types_check(equations_set_setup%FIELD,[field_u_variable_type],err,error,*999)
1329 CALL field_dimension_check(equations_set_setup%FIELD,field_u_variable_type,field_vector_dimension_type, &
1331 CALL field_data_type_check(equations_set_setup%FIELD,field_u_variable_type,field_dp_type,err,error,*999)
1332 CALL field_number_of_components_get(equations_set%GEOMETRY%GEOMETRIC_FIELD, field_u_variable_type, &
1333 & number_of_dimensions, err, error, *999)
1334 number_of_source_components = number_of_dimensions + 1
1335 CALL field_number_of_components_check(equations_set_setup%FIELD,field_u_variable_type, &
1336 & number_of_source_components,err,error,*999)
1339 CALL flagerror(
"Equations set source is not associated.",err,error,*999)
1342 equations_source=>equations_set%SOURCE
1343 IF(
ASSOCIATED(equations_source))
THEN 1344 IF(equations_source%SOURCE_FIELD_AUTO_CREATED)
THEN 1346 CALL field_create_finish(equations_source%SOURCE_FIELD,err,error,*999)
1348 CALL field_number_of_components_get(equations_set%GEOMETRY%GEOMETRIC_FIELD,field_u_variable_type, &
1349 & number_of_dimensions,err,error,*999)
1352 number_of_source_components = number_of_dimensions + 1
1354 number_of_source_components=0
1357 DO component_idx=1,number_of_source_components
1358 CALL field_component_values_initialise(equations_source%SOURCE_FIELD,field_u_variable_type, &
1359 & field_values_set_type,component_idx,0.0_dp,err,error,*999)
1363 CALL flagerror(
"Equations set source is not associated.",err,error,*999)
1366 local_error=
"The action type of "//
trim(
number_to_vstring(equations_set_setup%ACTION_TYPE,
"*",err,error))// &
1367 &
" for a setup type of "//
trim(
number_to_vstring(equations_set_setup%SETUP_TYPE,
"*",err,error))// &
1368 &
" is invalid for a standard, quasistatic or ALE Darcy equation." 1369 CALL flagerror(local_error,err,error,*999)
1377 SELECT CASE(equations_set%SPECIFICATION(3))
1382 SELECT CASE(equations_set_setup%ACTION_TYPE)
1384 equations_materials=>equations_set%MATERIALS
1385 IF(
ASSOCIATED(equations_materials) )
THEN 1386 IF( equations_materials%MATERIALS_FINISHED )
THEN 1391 CALL flagerror(
"Equations set materials has not been finished.",err,error,*999)
1394 CALL flagerror(
"Equations set materials is not associated.",err,error,*999)
1397 SELECT CASE(equations_set%SOLUTION_METHOD)
1399 SELECT CASE(equations_set%SPECIFICATION(3))
1407 equations_set_field_field=>equations_set%EQUATIONS_SET_FIELD%EQUATIONS_SET_FIELD_FIELD
1408 CALL field_parameter_set_data_get(equations_set_field_field,field_u_variable_type, &
1409 & field_values_set_type,equations_set_field_data,err,error,*999)
1410 imy_matrix = equations_set_field_data(1)
1411 ncompartments = equations_set_field_data(2)
1412 ALLOCATE(variable_types(2*ncompartments))
1413 DO num_var=1,ncompartments
1414 variable_types(2*num_var-1)=field_u_variable_type+(field_number_of_variable_subtypes*(num_var-1))
1415 variable_types(2*num_var)=field_deludeln_variable_type+(field_number_of_variable_subtypes*(num_var-1))
1423 SELECT CASE(equations%SPARSITY_TYPE)
1433 local_error=
"The equations matrices sparsity type of "// &
1435 CALL flagerror(local_error,err,error,*999)
1451 SELECT CASE(equations%SPARSITY_TYPE)
1461 local_error=
"The equations matrices sparsity type of "// &
1463 CALL flagerror(local_error,err,error,*999)
1468 CALL flagerror(
"Not implemented.",err,error,*999)
1470 CALL flagerror(
"Not implemented.",err,error,*999)
1472 CALL flagerror(
"Not implemented.",err,error,*999)
1474 CALL flagerror(
"Not implemented.",err,error,*999)
1476 CALL flagerror(
"Not implemented.",err,error,*999)
1478 local_error=
"The solution method of "//
trim(
number_to_vstring(equations_set%SOLUTION_METHOD,
"*",err,error))// &
1480 CALL flagerror(local_error,err,error,*999)
1484 local_error=
"The action type of "//
trim(
number_to_vstring(equations_set_setup%ACTION_TYPE,
"*",err,error))// &
1485 &
" for a setup type of "//
trim(
number_to_vstring(equations_set_setup%SETUP_TYPE,
"*",err,error))// &
1486 &
" is invalid for a standard Darcy equation." 1487 CALL flagerror(local_error,err,error,*999)
1494 SELECT CASE(equations_set_setup%ACTION_TYPE)
1496 equations_materials=>equations_set%MATERIALS
1497 IF(
ASSOCIATED(equations_materials) )
THEN 1498 IF( equations_materials%MATERIALS_FINISHED )
THEN 1503 CALL flagerror(
"Equations set materials has not been finished.",err,error,*999)
1506 CALL flagerror(
"Equations set materials is not associated.",err,error,*999)
1509 SELECT CASE(equations_set%SOLUTION_METHOD)
1517 SELECT CASE(equations_set%SPECIFICATION(3))
1530 SELECT CASE(equations%SPARSITY_TYPE)
1540 local_error=
"The equations matrices sparsity type of "// &
1542 CALL flagerror(local_error,err,error,*999)
1546 CALL flagerror(
"Not implemented.",err,error,*999)
1548 CALL flagerror(
"Not implemented.",err,error,*999)
1550 CALL flagerror(
"Not implemented.",err,error,*999)
1552 CALL flagerror(
"Not implemented.",err,error,*999)
1554 CALL flagerror(
"Not implemented.",err,error,*999)
1556 local_error=
"The solution method of "//
trim(
number_to_vstring(equations_set%SOLUTION_METHOD,
"*",err,error))// &
1558 CALL flagerror(local_error,err,error,*999)
1561 local_error=
"The action type of "//
trim(
number_to_vstring(equations_set_setup%ACTION_TYPE,
"*",err,error))// &
1562 &
" for a setup type of "//
trim(
number_to_vstring(equations_set_setup%SETUP_TYPE,
"*",err,error))// &
1563 &
" is invalid for a quasistatic Darcy equation." 1564 CALL flagerror(local_error,err,error,*999)
1572 SELECT CASE(equations_set_setup%ACTION_TYPE)
1574 equations_materials=>equations_set%MATERIALS
1575 IF(
ASSOCIATED(equations_materials))
THEN 1576 IF(equations_materials%MATERIALS_FINISHED)
THEN 1581 CALL flagerror(
"Equations set materials has not been finished.",err,error,*999)
1584 CALL flagerror(
"Equations materials is not associated.",err,error,*999)
1587 SELECT CASE(equations_set%SOLUTION_METHOD)
1599 SELECT CASE(equations_set%SPECIFICATION(3))
1609 equations_set_field_field=>equations_set%EQUATIONS_SET_FIELD%EQUATIONS_SET_FIELD_FIELD
1610 CALL field_parameter_set_data_get(equations_set_field_field,field_u_variable_type, &
1611 & field_values_set_type,equations_set_field_data,err,error,*999)
1612 imy_matrix = equations_set_field_data(1)
1613 ncompartments = equations_set_field_data(2)
1615 ALLOCATE(variable_types(2*ncompartments+2))
1616 ALLOCATE(variable_u_types(ncompartments-1))
1617 DO num_var=1,ncompartments+1
1618 variable_types(2*num_var-1)=field_u_variable_type+(field_number_of_variable_subtypes*(num_var-1))
1619 variable_types(2*num_var)=field_deludeln_variable_type+(field_number_of_variable_subtypes*(num_var-1))
1622 DO num_var=2,ncompartments+1
1623 IF((num_var-1)/=imy_matrix)
THEN 1624 num_var_count=num_var_count+1
1625 variable_u_types(num_var_count)=variable_types(2*num_var-1)
1652 SELECT CASE(equations%SPARSITY_TYPE)
1663 ALLOCATE(coupling_matrix_storage_type(ncompartments-1))
1664 ALLOCATE(coupling_matrix_structure_type(ncompartments-1))
1665 DO num_var=1,ncompartments-1
1670 & coupling_matrix_storage_type, &
1673 coupling_matrix_structure_type,err,error,*999)
1676 local_error=
"The equations matrices sparsity type of "// &
1678 CALL flagerror(local_error,err,error,*999)
1684 & err,error))//
" is invalid." 1685 CALL flagerror(local_error,err,error,*999)
1688 local_error=
"The action type of "//
trim(
number_to_vstring(equations_set_setup%ACTION_TYPE,
"*",err,error))// &
1689 &
" for a setup type of "//
trim(
number_to_vstring(equations_set_setup%SETUP_TYPE,
"*",err,error))// &
1690 &
" is invalid for a Darcy equation." 1691 CALL flagerror(local_error,err,error,*999)
1697 local_error=
"The equation set subtype of "//
trim(
number_to_vstring(equations_set%SPECIFICATION(3),
"*",err,error))// &
1699 &
" is invalid for a Darcy equation." 1700 CALL flagerror(local_error,err,error,*999)
1707 local_error=
"The setup type of "//
trim(
number_to_vstring(equations_set_setup%SETUP_TYPE,
"*",err,error))// &
1708 &
" is invalid for a standard, quasistatic, ALE or dynamic Darcy equation." 1709 CALL flagerror(local_error,err,error,*999)
1713 local_error=
"The equations set subtype of "//
trim(
number_to_vstring(equations_set%SPECIFICATION(3),
"*",err,error))// &
1714 &
" does not equal a standard, quasistatic, ALE or dynamic Darcy equation subtype." 1715 CALL flagerror(local_error,err,error,*999)
1718 CALL flagerror(
"Equations set is not associated.",err,error,*999)
1721 exits(
"DARCY_EQUATION_EQUATIONS_SET_SETUP")
1723 999 errorsexits(
"DARCY_EQUATION_EQUATIONS_SET_SETUP",err,error)
1736 INTEGER(INTG),
INTENT(IN) :: ELEMENT_NUMBER
1737 INTEGER(INTG),
INTENT(OUT) :: ERR
1741 INTEGER(INTG) :: FIELD_VAR_TYPE,ng,mh,mhs,mi,ms,nh,nhs,ni,ns,idxdim,num_var_count,idx_tensor
1742 INTEGER(INTG) :: my_compartment,Ncompartments,imatrix
1743 INTEGER(INTG) :: component_idx,xi_idx,derivative_idx
1744 INTEGER(INTG) :: MESH_COMPONENT_NUMBER, global_element_idx
1745 INTEGER(INTG) :: MESH_COMPONENT_1, MESH_COMPONENT_2
1746 INTEGER(INTG) :: NDOFS, NUMBER_OF_VEL_PRESS_COMPONENTS
1747 INTEGER(INTG) :: FIELD_VAR_TYPES(99)
1748 INTEGER(INTG) :: EQUATIONS_SET_SUBTYPE
1749 INTEGER(INTG),
POINTER :: EQUATIONS_SET_FIELD_DATA(:)
1750 REAL(DP) :: RWG,SUM,PGMSI(3),PGNSI(3),PGM,PGN,COUPLING_PARAM
1751 TYPE(
basis_type),
POINTER :: DEPENDENT_BASIS,GEOMETRIC_BASIS
1752 TYPE(
basis_type),
POINTER :: DEPENDENT_BASIS_1, DEPENDENT_BASIS_2
1762 TYPE(
field_type),
POINTER :: DEPENDENT_FIELD,GEOMETRIC_FIELD,MATERIALS_FIELD,EQUATIONS_SET_FIELD
1772 REAL(DP),
ALLOCATABLE :: PRESSURE_COEFF(:),PRESSURE(:),GRAD_PRESSURE(:,:)
1782 REAL(DP):: SOURCE,INTER_COMP_SOURCE,INTER_COMP_PERM_1,INTER_COMP_PERM_2
1783 REAL(DP):: BETA_PARAM, P_SINK_PARAM
1785 REAL(DP):: PERM_OVER_VIS_PARAM, DARCY_RHO_0_F
1786 REAL(DP):: PERM_TENSOR_OVER_VIS(3,3), VIS_OVER_PERM_TENSOR(3,3), Jmat
1787 REAL(DP):: X(3), ARG(3), L, FACT
1788 REAL(DP):: LM_PRESSURE,GRAD_LM_PRESSURE(3)
1790 REAL(DP):: DXDY(3,3), DXDXI(3,3), DYDXI(3,3), DXIDY(3,3)
1791 REAL(DP):: Jxy, Jyxi
1793 REAL(DP):: Mfact, bfact, p0fact
1800 LOGICAL :: STABILIZED
1805 darcy%LENGTH = 10.0_dp
1810 darcy%ANALYTIC = .false.
1812 enters(
"DARCY_EQUATION_FINITE_ELEMENT_CALCULATE",err,error,*999)
1817 NULLIFY(dependent_basis,geometric_basis)
1818 NULLIFY(dependent_basis_1, dependent_basis_2)
1820 NULLIFY(equations_mapping)
1821 NULLIFY(linear_mapping)
1822 NULLIFY(dynamic_mapping)
1823 NULLIFY(equations_matrices)
1824 NULLIFY(linear_matrices)
1825 NULLIFY(dynamic_matrices)
1827 NULLIFY(stiffness_matrix, damping_matrix)
1828 NULLIFY(dependent_field,geometric_field,materials_field)
1829 NULLIFY(field_variable)
1830 NULLIFY(quadrature_scheme)
1831 NULLIFY(quadrature_scheme_1, quadrature_scheme_2)
1832 NULLIFY(geometric_interpolated_point,materials_interpolated_point)
1833 NULLIFY(reference_geometric_interpolated_point)
1834 NULLIFY(elasticity_dependent_interpolation_parameters)
1835 NULLIFY(elasticity_dependent_interpolated_point)
1836 NULLIFY(decomposition,mesh_element)
1837 NULLIFY(boundary_conditions,boundary_conditions_variable)
1838 NULLIFY(source_vector,source_field)
1839 NULLIFY(equations_set_field_data)
1841 IF(
ASSOCIATED(equations_set))
THEN 1842 equations=>equations_set%EQUATIONS
1843 IF(
ASSOCIATED(equations))
THEN 1844 IF(.NOT.
ALLOCATED(equations_set%SPECIFICATION))
THEN 1845 CALL flagerror(
"Equations set specification is not allocated.",err,error,*999)
1846 ELSE IF(
SIZE(equations_set%SPECIFICATION,1)/=3)
THEN 1847 CALL flagerror(
"Equations set specification must have three entries for a Darcy type equations set.", &
1850 equations_set_subtype=equations_set%SPECIFICATION(3)
1851 SELECT CASE(equations_set_subtype)
1859 dependent_field=>equations%INTERPOLATION%DEPENDENT_FIELD
1860 geometric_field=>equations%INTERPOLATION%GEOMETRIC_FIELD
1861 materials_field=>equations%INTERPOLATION%MATERIALS_FIELD
1864 source_field=>equations%INTERPOLATION%SOURCE_FIELD
1867 equations_matrices=>equations%EQUATIONS_MATRICES
1868 rhs_vector=>equations_matrices%RHS_VECTOR
1869 equations_mapping=>equations%EQUATIONS_MAPPING
1873 source_vector=>equations_matrices%SOURCE_VECTOR
1874 source_vector%ELEMENT_VECTOR%VECTOR = 0.0_dp
1877 SELECT CASE(equations_set_subtype)
1880 linear_matrices=>equations_matrices%LINEAR_MATRICES
1881 stiffness_matrix=>linear_matrices%MATRICES(1)%PTR
1883 linear_mapping=>equations_mapping%LINEAR_MAPPING
1884 field_variable=>linear_mapping%EQUATIONS_MATRIX_TO_VAR_MAPS(1)%VARIABLE
1885 field_var_type=field_variable%VARIABLE_TYPE
1887 stiffness_matrix%ELEMENT_MATRIX%MATRIX=0.0_dp
1892 dynamic_matrices=>equations_matrices%DYNAMIC_MATRICES
1893 stiffness_matrix=>dynamic_matrices%MATRICES(1)%PTR
1894 damping_matrix=>dynamic_matrices%MATRICES(2)%PTR
1896 dynamic_mapping=>equations_mapping%DYNAMIC_MAPPING
1897 field_variable=>dynamic_mapping%EQUATIONS_MATRIX_TO_VAR_MAPS(1)%VARIABLE
1898 field_var_type=field_variable%VARIABLE_TYPE
1900 stiffness_matrix%ELEMENT_MATRIX%MATRIX=0.0_dp
1901 damping_matrix%ELEMENT_MATRIX%MATRIX=0.0_dp
1904 decomposition => dependent_field%DECOMPOSITION
1905 mesh_component_number = decomposition%MESH_COMPONENT_NUMBER
1906 global_element_idx = decomposition%DOMAIN(mesh_component_number)%PTR%MAPPINGS%ELEMENTS% &
1907 & local_to_global_map(element_number)
1908 mesh_element => decomposition%MESH%TOPOLOGY(mesh_component_number)%PTR%ELEMENTS%ELEMENTS(global_element_idx)
1911 equations_set_field=>equations_set%EQUATIONS_SET_FIELD%EQUATIONS_SET_FIELD_FIELD
1912 CALL field_parameter_set_data_get(equations_set_field,field_u_variable_type, &
1913 & field_values_set_type,equations_set_field_data,err,error,*999)
1915 my_compartment = equations_set_field_data(1)
1916 ncompartments = equations_set_field_data(2)
1917 linear_matrices=>equations_matrices%LINEAR_MATRICES
1918 stiffness_matrix=>linear_matrices%MATRICES(1)%PTR
1920 linear_mapping=>equations_mapping%LINEAR_MAPPING
1921 field_variable=>linear_mapping%EQUATIONS_MATRIX_TO_VAR_MAPS(1)%VARIABLE
1922 field_var_type=field_variable%VARIABLE_TYPE
1924 stiffness_matrix%ELEMENT_MATRIX%MATRIX=0.0_dp
1928 equations_set_field=>equations_set%EQUATIONS_SET_FIELD%EQUATIONS_SET_FIELD_FIELD
1929 CALL field_parameter_set_data_get(equations_set_field,field_u_variable_type, &
1930 & field_values_set_type,equations_set_field_data,err,error,*999)
1932 my_compartment = equations_set_field_data(1)
1933 ncompartments = equations_set_field_data(2)
1937 linear_matrices=>equations_matrices%LINEAR_MATRICES
1938 linear_mapping=>equations_mapping%LINEAR_MAPPING
1947 dynamic_matrices=>equations_matrices%DYNAMIC_MATRICES
1948 stiffness_matrix=>dynamic_matrices%MATRICES(1)%PTR
1949 damping_matrix=>dynamic_matrices%MATRICES(2)%PTR
1953 dynamic_mapping=>equations_mapping%DYNAMIC_MAPPING
1954 field_variable=>dynamic_mapping%EQUATIONS_MATRIX_TO_VAR_MAPS(1)%VARIABLE
1955 field_var_type=field_variable%VARIABLE_TYPE
1957 stiffness_matrix%ELEMENT_MATRIX%MATRIX=0.0_dp
1958 damping_matrix%ELEMENT_MATRIX%MATRIX=0.0_dp
1960 equations_set_field=>equations_set%EQUATIONS_SET_FIELD%EQUATIONS_SET_FIELD_FIELD
1961 CALL field_parameter_set_data_get(equations_set_field,field_u_variable_type, &
1962 & field_values_set_type,equations_set_field_data,err,error,*999)
1964 my_compartment = equations_set_field_data(1)
1965 ncompartments = equations_set_field_data(2)
1968 linear_matrices=>equations_matrices%LINEAR_MATRICES
1969 linear_mapping=>equations_mapping%LINEAR_MAPPING
1972 DO imatrix = 1,ncompartments
1973 IF(imatrix/=my_compartment)
THEN 1974 num_var_count=num_var_count+1
1975 coupling_matrices(num_var_count)%PTR=>linear_matrices%MATRICES(num_var_count)%PTR
1976 field_variables(num_var_count)%PTR=>linear_mapping%EQUATIONS_MATRIX_TO_VAR_MAPS(num_var_count)%VARIABLE
1977 field_var_types(num_var_count)=field_variables(num_var_count)%PTR%VARIABLE_TYPE
1978 coupling_matrices(num_var_count)%PTR%ELEMENT_MATRIX%MATRIX=0.0_dp
1982 dynamic_matrices=>equations_matrices%DYNAMIC_MATRICES
1983 stiffness_matrix=>dynamic_matrices%MATRICES(1)%PTR
1984 damping_matrix=>dynamic_matrices%MATRICES(2)%PTR
1986 dynamic_mapping=>equations_mapping%DYNAMIC_MAPPING
1987 field_variable=>dynamic_mapping%EQUATIONS_MATRIX_TO_VAR_MAPS(1)%VARIABLE
1988 field_var_type=field_variable%VARIABLE_TYPE
1990 stiffness_matrix%ELEMENT_MATRIX%MATRIX=0.0_dp
1991 damping_matrix%ELEMENT_MATRIX%MATRIX=0.0_dp
1993 ALLOCATE(pressure_coeff(ncompartments))
1995 ALLOCATE(pressure(ncompartments))
1996 ALLOCATE(grad_pressure(3,ncompartments))
1998 grad_pressure = 0.0_dp
1999 pressure_coeff(1)=0.25_dp
2000 pressure_coeff(2)=0.25_dp
2001 pressure_coeff(3)=0.25_dp
2002 pressure_coeff(4)=0.25_dp
2007 geometric_basis=>geometric_field%DECOMPOSITION%DOMAIN(geometric_field%DECOMPOSITION%MESH_COMPONENT_NUMBER)%PTR% &
2008 & topology%ELEMENTS%ELEMENTS(element_number)%BASIS
2009 dependent_basis=>dependent_field%DECOMPOSITION%DOMAIN(dependent_field%DECOMPOSITION%MESH_COMPONENT_NUMBER)%PTR% &
2010 & topology%ELEMENTS%ELEMENTS(element_number)%BASIS
2014 CALL field_interpolation_parameters_element_get(field_values_set_type,element_number,equations%INTERPOLATION% &
2015 & geometric_interp_parameters(field_u_variable_type)%PTR,err,error,*999)
2016 CALL field_interpolation_parameters_element_get(field_values_set_type,element_number,equations%INTERPOLATION% &
2017 & materials_interp_parameters(field_u_variable_type)%PTR,err,error,*999)
2019 CALL field_interpolation_parameters_element_get(field_values_set_type,element_number,equations%INTERPOLATION% &
2020 & materials_interp_parameters(field_v_variable_type)%PTR,err,error,*999)
2021 CALL field_interpolation_parameters_element_get(field_values_set_type,element_number,equations%INTERPOLATION% &
2022 & materials_interp_parameters(field_u1_variable_type)%PTR,err,error,*999)
2027 CALL field_interpolation_parameters_element_get(field_values_set_type,element_number,equations%INTERPOLATION% &
2028 & source_interp_parameters(field_u_variable_type)%PTR,err,error,*999)
2033 elasticity_dependent_interpolation_parameters=>equations%INTERPOLATION% &
2034 & dependent_interp_parameters(field_u_variable_type)%PTR
2035 CALL field_interpolation_parameters_element_get(field_values_set_type,element_number, &
2036 & elasticity_dependent_interpolation_parameters,err,error,*999)
2037 elasticity_dependent_interpolated_point=>equations%INTERPOLATION% &
2038 & dependent_interp_point(field_u_variable_type)%PTR
2042 SELECT CASE(equations_set_subtype)
2045 number_of_vel_press_components = field_variable%NUMBER_OF_COMPONENTS - 1
2047 number_of_vel_press_components = field_variable%NUMBER_OF_COMPONENTS
2054 SELECT CASE(equations_set_subtype)
2058 IF( mesh_element%BOUNDARY_ELEMENT )
THEN 2066 DO ng=1,quadrature_scheme%NUMBER_OF_GAUSS
2076 CALL field_interpolation_parameters_element_get(field_initial_values_set_type,element_number, &
2077 & equations%INTERPOLATION%GEOMETRIC_INTERP_PARAMETERS(field_u_variable_type)%PTR,err,error,*999)
2078 reference_geometric_interpolated_point => equations%INTERPOLATION%GEOMETRIC_INTERP_POINT(field_u_variable_type)%PTR
2080 & reference_geometric_interpolated_point,err,error,*999)
2082 DO component_idx=1,dependent_basis%NUMBER_OF_XI
2083 DO xi_idx=1,dependent_basis%NUMBER_OF_XI
2085 dydxi(component_idx,xi_idx)=reference_geometric_interpolated_point%VALUES(component_idx,derivative_idx)
2090 CALL field_interpolation_parameters_element_get(field_values_set_type,element_number, &
2091 & equations%INTERPOLATION%GEOMETRIC_INTERP_PARAMETERS(field_u_variable_type)%PTR,err,error,*999)
2092 geometric_interpolated_point => equations%INTERPOLATION%GEOMETRIC_INTERP_POINT(field_u_variable_type)%PTR
2094 & geometric_interpolated_point,err,error,*999)
2095 CALL field_interpolated_point_metrics_calculate(geometric_basis%NUMBER_OF_XI, &
2096 & equations%INTERPOLATION%GEOMETRIC_INTERP_POINT_METRICS(field_u_variable_type)%PTR,err,error,*999)
2098 DO component_idx=1,dependent_basis%NUMBER_OF_XI
2099 DO xi_idx=1,dependent_basis%NUMBER_OF_XI
2101 dxdxi(component_idx,xi_idx)=geometric_interpolated_point%VALUES(component_idx,derivative_idx)
2106 CALL invert(dydxi,dxidy,jyxi,err,error,*999)
2110 IF( abs(jxy) < 1.0e-10_dp )
THEN 2111 local_error=
"DARCY_EQUATION_FINITE_ELEMENT_CALCULATE: Jacobian Jxy is smaller than 1.0E-10_DP." 2112 CALL flagerror(local_error,err,error,*999)
2123 geometric_interpolated_point => equations%INTERPOLATION%GEOMETRIC_INTERP_POINT(field_u_variable_type)%PTR
2125 & geometric_interpolated_point,err,error,*999)
2126 CALL field_interpolated_point_metrics_calculate(geometric_basis%NUMBER_OF_XI, &
2127 & equations%INTERPOLATION%GEOMETRIC_INTERP_POINT_METRICS(field_u_variable_type)%PTR,err,error,*999)
2130 CALL field_interpolation_parameters_element_get(field_values_set_type,element_number, &
2131 & equations%INTERPOLATION%GEOMETRIC_INTERP_PARAMETERS(field_u_variable_type)%PTR,err,error,*999)
2133 & geometric_interpolated_point,err,error,*999)
2142 materials_interpolated_point => equations%INTERPOLATION%MATERIALS_INTERP_POINT(field_u_variable_type)%PTR
2144 & materials_interpolated_point,err,error,*999)
2148 & materials_interp_point(field_v_variable_type)%PTR,err,error,*999)
2151 & materials_interp_point(field_u1_variable_type)%PTR,err,error,*999)
2157 & source_interp_point(field_u_variable_type)%PTR,err,error,*999)
2160 SELECT CASE(equations_set_subtype)
2163 perm_tensor_over_vis=0.0_dp
2164 perm_tensor_over_vis(1,1) = materials_interpolated_point%VALUES(2,
no_part_deriv)
2165 perm_tensor_over_vis(2,2) = materials_interpolated_point%VALUES(2,
no_part_deriv)
2166 perm_tensor_over_vis(3,3) = materials_interpolated_point%VALUES(2,
no_part_deriv)
2168 perm_tensor_over_vis=perm_tensor_over_vis*materials_interpolated_point%VALUES(1,
no_part_deriv)
2171 perm_tensor_over_vis(1,1) = materials_interpolated_point%VALUES(2,
no_part_deriv)
2172 perm_tensor_over_vis(1,2) = materials_interpolated_point%VALUES(3,
no_part_deriv)
2173 perm_tensor_over_vis(1,3) = materials_interpolated_point%VALUES(4,
no_part_deriv)
2174 perm_tensor_over_vis(2,2) = materials_interpolated_point%VALUES(5,
no_part_deriv)
2175 perm_tensor_over_vis(2,3) = materials_interpolated_point%VALUES(6,
no_part_deriv)
2176 perm_tensor_over_vis(3,3) = materials_interpolated_point%VALUES(7,
no_part_deriv)
2178 perm_tensor_over_vis(2,1) = perm_tensor_over_vis(1,2)
2179 perm_tensor_over_vis(3,1) = perm_tensor_over_vis(1,3)
2180 perm_tensor_over_vis(3,2) = perm_tensor_over_vis(2,3)
2186 & materials_interpolated_point%VALUES(1,
no_part_deriv),err,error,*999)
2188 & materials_interpolated_point%VALUES(2,
no_part_deriv),err,error,*999)
2194 jmat =
determinant(perm_tensor_over_vis,err,error)
2196 CALL invert(perm_tensor_over_vis,vis_over_perm_tensor,jmat,err,error,*999)
2198 vis_over_perm_tensor = 0.0_dp
2200 vis_over_perm_tensor(idx_tensor,idx_tensor) = 1.0e10_dp
2208 beta_param = -
darcy%PERM_OVER_VIS * (2.0_dp *
pi /
darcy%LENGTH) * (2.0_dp *
pi /
darcy%LENGTH)
2209 p_sink_param =
darcy%P_SINK
2214 & elasticity_dependent_interpolated_point,err,error,*999)
2218 lm_pressure = -elasticity_dependent_interpolated_point%VALUES(4,
no_part_deriv)
2219 DO xi_idx=1,dependent_basis%NUMBER_OF_XI
2222 grad_lm_pressure(xi_idx) = -elasticity_dependent_interpolated_point%VALUES(4,derivative_idx)
2228 write(*,*)
'NEED CONSTITUTIVE LAWS HERE!!!! THE FOLLOWING IS PLACEHOLDER ONLY!' 2231 & elasticity_dependent_interpolated_point,err,error,*999)
2233 lm_pressure = -elasticity_dependent_interpolated_point%VALUES(4,
no_part_deriv)
2234 DO xi_idx=1,dependent_basis%NUMBER_OF_XI
2237 grad_lm_pressure(xi_idx) = -elasticity_dependent_interpolated_point%VALUES(4,derivative_idx)
2246 DO imatrix=1,ncompartments
2248 pressure(imatrix) = pressure_coeff(imatrix)*lm_pressure
2249 DO xi_idx=1,dependent_basis%NUMBER_OF_XI
2250 grad_pressure(xi_idx,imatrix) = pressure_coeff(imatrix)*grad_lm_pressure(xi_idx)
2260 DO mh=1,field_variable%NUMBER_OF_COMPONENTS
2262 mesh_component_1 = field_variable%COMPONENTS(mh)%MESH_COMPONENT_NUMBER
2263 dependent_basis_1 => dependent_field%DECOMPOSITION%DOMAIN(mesh_component_1)%PTR% &
2264 & topology%ELEMENTS%ELEMENTS(element_number)%BASIS
2265 quadrature_scheme_1 => dependent_basis_1%QUADRATURE% &
2267 rwg = equations%INTERPOLATION%GEOMETRIC_INTERP_POINT_METRICS(field_u_variable_type)%PTR%JACOBIAN * &
2268 & quadrature_scheme_1%GAUSS_WEIGHTS(ng)
2270 DO ms=1,dependent_basis_1%NUMBER_OF_ELEMENT_PARAMETERS
2275 IF(stiffness_matrix%UPDATE_MATRIX)
THEN 2279 DO nh=1,field_variable%NUMBER_OF_COMPONENTS
2281 mesh_component_2 = field_variable%COMPONENTS(nh)%MESH_COMPONENT_NUMBER
2282 dependent_basis_2 => dependent_field%DECOMPOSITION%DOMAIN(mesh_component_2)%PTR% &
2283 & topology%ELEMENTS%ELEMENTS(element_number)%BASIS
2285 quadrature_scheme_2 => dependent_basis_2%QUADRATURE% &
2290 DO ns=1,dependent_basis_2%NUMBER_OF_ELEMENT_PARAMETERS
2293 SELECT CASE(equations_set_subtype)
2299 IF(mh==nh.AND.nh<field_variable%NUMBER_OF_COMPONENTS)
THEN 2303 pgm=quadrature_scheme_1%GAUSS_BASIS_FNS(ms,
no_part_deriv,ng)
2304 pgn=quadrature_scheme_2%GAUSS_BASIS_FNS(ns,
no_part_deriv,ng)
2306 sum = sum + vis_over_perm_tensor( mh, nh ) * pgm * pgn
2310 stiffness_matrix%ELEMENT_MATRIX%MATRIX(mhs,nhs) = stiffness_matrix%ELEMENT_MATRIX%MATRIX(mhs,nhs) + &
2314 ELSE IF(mh==field_variable%NUMBER_OF_COMPONENTS.AND.nh<field_variable%NUMBER_OF_COMPONENTS)
THEN 2318 pgm=quadrature_scheme_1%GAUSS_BASIS_FNS(ms,
no_part_deriv,ng)
2319 pgn=quadrature_scheme_2%GAUSS_BASIS_FNS(ns,
no_part_deriv,ng)
2321 DO mi=1,dependent_basis_1%NUMBER_OF_XI
2325 DO ni=1,dependent_basis_2%NUMBER_OF_XI
2329 DO ni=1,dependent_basis_2%NUMBER_OF_XI
2330 sum = sum + pgm * pgnsi(ni) * &
2331 & equations%INTERPOLATION%GEOMETRIC_INTERP_POINT_METRICS(field_u_variable_type)%PTR%DXI_DX(ni,nh)
2334 stiffness_matrix%ELEMENT_MATRIX%MATRIX(mhs,nhs) = stiffness_matrix%ELEMENT_MATRIX%MATRIX(mhs,nhs) + &
2339 IF(damping_matrix%UPDATE_MATRIX)
THEN 2341 IF(mh==nh.AND.nh==field_variable%NUMBER_OF_COMPONENTS)
THEN 2342 pgm=quadrature_scheme_1%GAUSS_BASIS_FNS(ms,
no_part_deriv,ng)
2343 pgn=quadrature_scheme_2%GAUSS_BASIS_FNS(ns,
no_part_deriv,ng)
2348 sum = pgm * pgn / (jxy * darcy_rho_0_f)
2350 damping_matrix%ELEMENT_MATRIX%MATRIX(mhs,nhs) = damping_matrix%ELEMENT_MATRIX%MATRIX(mhs,nhs) + &
2372 IF(mh==nh.AND.nh<field_variable%NUMBER_OF_COMPONENTS)
THEN 2376 pgm=quadrature_scheme_1%GAUSS_BASIS_FNS(ms,
no_part_deriv,ng)
2377 pgn=quadrature_scheme_2%GAUSS_BASIS_FNS(ns,
no_part_deriv,ng)
2379 sum = sum + vis_over_perm_tensor( mh, nh ) * pgm * pgn
2383 stiffness_matrix%ELEMENT_MATRIX%MATRIX(mhs,nhs) = stiffness_matrix%ELEMENT_MATRIX%MATRIX(mhs,nhs) + &
2387 ELSE IF(mh==field_variable%NUMBER_OF_COMPONENTS.AND.nh<field_variable%NUMBER_OF_COMPONENTS)
THEN 2391 pgm=quadrature_scheme_1%GAUSS_BASIS_FNS(ms,
no_part_deriv,ng)
2392 pgn=quadrature_scheme_2%GAUSS_BASIS_FNS(ns,
no_part_deriv,ng)
2394 DO mi=1,dependent_basis_1%NUMBER_OF_XI
2398 DO ni=1,dependent_basis_2%NUMBER_OF_XI
2402 DO ni=1,dependent_basis_2%NUMBER_OF_XI
2403 sum = sum + pgm * pgnsi(ni) * &
2404 & equations%INTERPOLATION%GEOMETRIC_INTERP_POINT_METRICS(field_u_variable_type)%PTR%DXI_DX(ni,nh)
2407 stiffness_matrix%ELEMENT_MATRIX%MATRIX(mhs,nhs) = stiffness_matrix%ELEMENT_MATRIX%MATRIX(mhs,nhs) + &
2412 IF(damping_matrix%UPDATE_MATRIX)
THEN 2414 IF(mh==nh.AND.nh==field_variable%NUMBER_OF_COMPONENTS)
THEN 2415 pgm=quadrature_scheme_1%GAUSS_BASIS_FNS(ms,
no_part_deriv,ng)
2416 pgn=quadrature_scheme_2%GAUSS_BASIS_FNS(ns,
no_part_deriv,ng)
2421 sum = pgm * pgn / (jxy * darcy_rho_0_f)
2423 damping_matrix%ELEMENT_MATRIX%MATRIX(mhs,nhs) = damping_matrix%ELEMENT_MATRIX%MATRIX(mhs,nhs) + &
2448 IF(mh==nh.AND.nh<number_of_vel_press_components)
THEN 2452 pgm=quadrature_scheme_1%GAUSS_BASIS_FNS(ms,
no_part_deriv,ng)
2453 pgn=quadrature_scheme_2%GAUSS_BASIS_FNS(ns,
no_part_deriv,ng)
2455 sum = sum + vis_over_perm_tensor( mh, nh ) * pgm * pgn
2459 IF( stabilized )
THEN 2460 sum = sum - 0.5_dp * vis_over_perm_tensor( mh, nh ) * pgm * pgn
2463 stiffness_matrix%ELEMENT_MATRIX%MATRIX(mhs,nhs) = stiffness_matrix%ELEMENT_MATRIX%MATRIX(mhs,nhs) + &
2468 ELSE IF(mh<number_of_vel_press_components.AND.nh==number_of_vel_press_components)
THEN 2472 pgm=quadrature_scheme_1%GAUSS_BASIS_FNS(ms,
no_part_deriv,ng)
2473 pgn=quadrature_scheme_2%GAUSS_BASIS_FNS(ns,
no_part_deriv,ng)
2475 DO mi=1,dependent_basis_1%NUMBER_OF_XI
2479 DO ni=1,dependent_basis_2%NUMBER_OF_XI
2483 DO mi=1,dependent_basis_1%NUMBER_OF_XI
2484 sum = sum - pgmsi(mi) * pgn * &
2485 & equations%INTERPOLATION%GEOMETRIC_INTERP_POINT_METRICS(field_u_variable_type)%PTR%DXI_DX(mi,mh)
2488 IF( stabilized )
THEN 2489 DO ni=1,dependent_basis_2%NUMBER_OF_XI
2490 sum = sum - 0.5_dp * pgm * pgnsi(ni) * &
2491 & equations%INTERPOLATION%GEOMETRIC_INTERP_POINT_METRICS(field_u_variable_type)%PTR%DXI_DX(ni,mh)
2495 stiffness_matrix%ELEMENT_MATRIX%MATRIX(mhs,nhs) = stiffness_matrix%ELEMENT_MATRIX%MATRIX(mhs,nhs) + &
2500 ELSE IF(mh==number_of_vel_press_components.AND.nh<number_of_vel_press_components)
THEN 2504 pgm=quadrature_scheme_1%GAUSS_BASIS_FNS(ms,
no_part_deriv,ng)
2505 pgn=quadrature_scheme_2%GAUSS_BASIS_FNS(ns,
no_part_deriv,ng)
2507 DO mi=1,dependent_basis_1%NUMBER_OF_XI
2511 DO ni=1,dependent_basis_2%NUMBER_OF_XI
2515 DO ni=1,dependent_basis_2%NUMBER_OF_XI
2516 sum = sum + pgm * pgnsi(ni) * &
2517 & equations%INTERPOLATION%GEOMETRIC_INTERP_POINT_METRICS(field_u_variable_type)%PTR%DXI_DX(ni,nh)
2520 IF( stabilized )
THEN 2521 DO mi=1,dependent_basis_1%NUMBER_OF_XI
2522 sum = sum + 0.5_dp * pgmsi(mi) * pgn * &
2523 & equations%INTERPOLATION%GEOMETRIC_INTERP_POINT_METRICS(field_u_variable_type)%PTR%DXI_DX(mi,nh)
2527 stiffness_matrix%ELEMENT_MATRIX%MATRIX(mhs,nhs) = stiffness_matrix%ELEMENT_MATRIX%MATRIX(mhs,nhs) + &
2532 ELSE IF(mh==nh.AND.nh==number_of_vel_press_components)
THEN 2536 DO mi=1,dependent_basis_1%NUMBER_OF_XI
2540 DO ni=1,dependent_basis_2%NUMBER_OF_XI
2544 IF( stabilized )
THEN 2545 DO idxdim =1,dependent_basis_1%NUMBER_OF_XI
2546 DO mi=1,dependent_basis_1%NUMBER_OF_XI
2547 DO ni=1,dependent_basis_2%NUMBER_OF_XI
2548 sum = sum + 0.5_dp * perm_tensor_over_vis( idxdim, idxdim ) * pgmsi(mi) * pgnsi(ni) * &
2549 & equations%INTERPOLATION%GEOMETRIC_INTERP_POINT_METRICS(field_u_variable_type)%PTR% &
2550 & dxi_dx(mi,idxdim) * equations%INTERPOLATION% &
2551 & geometric_interp_point_metrics(field_u_variable_type)%PTR%DXI_DX(ni,idxdim)
2557 IF(
darcy%TESTCASE == 3 )
THEN 2560 pgm=quadrature_scheme_1%GAUSS_BASIS_FNS(ms,
no_part_deriv,ng)
2561 pgn=quadrature_scheme_2%GAUSS_BASIS_FNS(ns,
no_part_deriv,ng)
2563 sum = sum + beta_param * pgm * pgn
2566 stiffness_matrix%ELEMENT_MATRIX%MATRIX(mhs,nhs) = stiffness_matrix%ELEMENT_MATRIX%MATRIX(mhs,nhs) + &
2572 & mh==field_variable%NUMBER_OF_COMPONENTS.AND.nh==number_of_vel_press_components)
THEN 2576 pgm=quadrature_scheme_1%GAUSS_BASIS_FNS(ms,
no_part_deriv,ng)
2577 pgn=quadrature_scheme_2%GAUSS_BASIS_FNS(ns,
no_part_deriv,ng)
2579 sum = sum - pgm * pgn / (mfact * ffact)
2581 stiffness_matrix%ELEMENT_MATRIX%MATRIX(mhs,nhs) = stiffness_matrix%ELEMENT_MATRIX%MATRIX(mhs,nhs) + &
2587 & mh==nh.AND.nh==field_variable%NUMBER_OF_COMPONENTS)
THEN 2591 pgm=quadrature_scheme_1%GAUSS_BASIS_FNS(ms,
no_part_deriv,ng)
2592 pgn=quadrature_scheme_2%GAUSS_BASIS_FNS(ns,
no_part_deriv,ng)
2594 sum = sum + pgm * pgn / darcy_rho_0_f
2596 stiffness_matrix%ELEMENT_MATRIX%MATRIX(mhs,nhs) = stiffness_matrix%ELEMENT_MATRIX%MATRIX(mhs,nhs) + &
2602 stiffness_matrix%ELEMENT_MATRIX%MATRIX(mhs,nhs) = 0.0_dp
2610 IF(damping_matrix%UPDATE_MATRIX)
THEN 2611 IF(mh==nh.AND.mh<number_of_vel_press_components)
THEN 2612 pgm=quadrature_scheme_1%GAUSS_BASIS_FNS(ms,
no_part_deriv,ng)
2613 pgn=quadrature_scheme_2%GAUSS_BASIS_FNS(ns,
no_part_deriv,ng)
2617 sum = pgm*pgn*darcy_rho_0_f
2619 damping_matrix%ELEMENT_MATRIX%MATRIX(mhs,nhs) = damping_matrix%ELEMENT_MATRIX%MATRIX(mhs,nhs) + &
2624 IF(damping_matrix%UPDATE_MATRIX)
THEN 2626 IF(mh==number_of_vel_press_components.AND.nh==field_variable%NUMBER_OF_COMPONENTS)
THEN 2627 pgm=quadrature_scheme_1%GAUSS_BASIS_FNS(ms,
no_part_deriv,ng)
2628 pgn=quadrature_scheme_2%GAUSS_BASIS_FNS(ns,
no_part_deriv,ng)
2632 sum = pgm * pgn / (jxy * darcy_rho_0_f)
2634 damping_matrix%ELEMENT_MATRIX%MATRIX(mhs,nhs) = damping_matrix%ELEMENT_MATRIX%MATRIX(mhs,nhs) + &
2649 IF(rhs_vector%UPDATE_VECTOR)
THEN 2651 SELECT CASE(equations_set_subtype)
2658 IF( mh<field_variable%NUMBER_OF_COMPONENTS )
THEN 2662 pgm = quadrature_scheme_1%GAUSS_BASIS_FNS(ms,
no_part_deriv,ng)
2665 DO mi=1,dependent_basis_1%NUMBER_OF_XI
2666 sum = sum - pgm * grad_lm_pressure(mi) * &
2667 & equations%INTERPOLATION%GEOMETRIC_INTERP_POINT_METRICS(field_u_variable_type)%PTR%DXI_DX(mi,mh)
2670 rhs_vector%ELEMENT_VECTOR%VECTOR(mhs) = rhs_vector%ELEMENT_VECTOR%VECTOR(mhs) + sum * rwg
2674 ELSE IF( mh==field_variable%NUMBER_OF_COMPONENTS )
THEN 2678 pgm=quadrature_scheme_1%GAUSS_BASIS_FNS(ms,
no_part_deriv,ng)
2683 sum = sum + pgm * source
2685 rhs_vector%ELEMENT_VECTOR%VECTOR(mhs) = rhs_vector%ELEMENT_VECTOR%VECTOR(mhs) + sum * rwg
2689 rhs_vector%ELEMENT_VECTOR%VECTOR(mhs) = 0.0_dp
2696 IF( mh<field_variable%NUMBER_OF_COMPONENTS )
THEN 2700 pgm = quadrature_scheme_1%GAUSS_BASIS_FNS(ms,
no_part_deriv,ng)
2704 DO mi=1,dependent_basis_1%NUMBER_OF_XI
2708 sum = sum - pgm * grad_pressure(mi,my_compartment) * &
2709 & equations%INTERPOLATION%GEOMETRIC_INTERP_POINT_METRICS(field_u_variable_type)%PTR%DXI_DX(mi,mh)
2712 rhs_vector%ELEMENT_VECTOR%VECTOR(mhs) = rhs_vector%ELEMENT_VECTOR%VECTOR(mhs) + sum * rwg
2716 ELSE IF( mh==field_variable%NUMBER_OF_COMPONENTS )
THEN 2720 pgm=quadrature_scheme_1%GAUSS_BASIS_FNS(ms,
no_part_deriv,ng)
2728 DO imatrix=1,ncompartments
2729 IF(imatrix/=my_compartment)
THEN 2731 inter_comp_perm_1=equations%INTERPOLATION%MATERIALS_INTERP_POINT(field_v_variable_type)%PTR% &
2733 inter_comp_perm_2=equations%INTERPOLATION%MATERIALS_INTERP_POINT(field_v_variable_type)%PTR% &
2736 inter_comp_source=-inter_comp_perm_1*pressure(my_compartment) + inter_comp_perm_2*pressure(imatrix)
2741 sum = sum + pgm * source + pgm * inter_comp_source
2743 rhs_vector%ELEMENT_VECTOR%VECTOR(mhs) = rhs_vector%ELEMENT_VECTOR%VECTOR(mhs) + sum * rwg
2747 rhs_vector%ELEMENT_VECTOR%VECTOR(mhs) = 0.0_dp
2755 IF( mh<number_of_vel_press_components )
THEN 2759 rhs_vector%ELEMENT_VECTOR%VECTOR(mhs) = rhs_vector%ELEMENT_VECTOR%VECTOR(mhs) + sum * rwg
2763 ELSE IF( mh==number_of_vel_press_components )
THEN 2767 pgm=quadrature_scheme_1%GAUSS_BASIS_FNS(ms,
no_part_deriv,ng)
2772 sum = sum + pgm * source
2774 rhs_vector%ELEMENT_VECTOR%VECTOR(mhs) = rhs_vector%ELEMENT_VECTOR%VECTOR(mhs) + sum * rwg
2779 & mh==field_variable%NUMBER_OF_COMPONENTS)
THEN 2783 pgm=quadrature_scheme_1%GAUSS_BASIS_FNS(ms,
no_part_deriv,ng)
2785 sum = sum - pgm * bfact * (1.0_dp - jxy)
2787 sum = sum - pgm * p0fact / (mfact * ffact)
2789 rhs_vector%ELEMENT_VECTOR%VECTOR(mhs) = rhs_vector%ELEMENT_VECTOR%VECTOR(mhs) + sum * rwg
2793 rhs_vector%ELEMENT_VECTOR%VECTOR(mhs) = 0.0_dp
2803 IF(
ASSOCIATED(source_vector))
THEN 2804 IF(source_vector%UPDATE_VECTOR)
THEN 2808 c_param=equations%INTERPOLATION%SOURCE_INTERP_POINT(field_u_variable_type)%PTR%VALUES(mh,
no_part_deriv)
2813 pgm=quadrature_scheme_1%GAUSS_BASIS_FNS(ms,
no_part_deriv,ng)
2814 sum = sum + pgm * c_param
2815 source_vector%ELEMENT_VECTOR%VECTOR(mhs) = source_vector%ELEMENT_VECTOR%VECTOR(mhs) + sum * rwg
2827 DO mh=1,field_variable%NUMBER_OF_COMPONENTS
2829 mesh_component_1 = field_variable%COMPONENTS(mh)%MESH_COMPONENT_NUMBER
2830 dependent_basis_1 => dependent_field%DECOMPOSITION%DOMAIN(mesh_component_1)%PTR% &
2831 & topology%ELEMENTS%ELEMENTS(element_number)%BASIS
2832 quadrature_scheme_1 => dependent_basis_1%QUADRATURE% &
2834 rwg = equations%INTERPOLATION%GEOMETRIC_INTERP_POINT_METRICS(field_u_variable_type)%PTR%JACOBIAN * &
2835 & quadrature_scheme_1%GAUSS_WEIGHTS(ng)
2837 DO ms=1,dependent_basis_1%NUMBER_OF_ELEMENT_PARAMETERS
2841 DO imatrix = 1,ncompartments
2842 IF(imatrix/=my_compartment)
THEN 2843 num_var_count=num_var_count+1
2847 IF(coupling_matrices(num_var_count)%PTR%UPDATE_MATRIX)
THEN 2851 DO nh=1,field_variables(num_var_count)%PTR%NUMBER_OF_COMPONENTS
2853 mesh_component_2 = field_variable%COMPONENTS(nh)%MESH_COMPONENT_NUMBER
2854 dependent_basis_2 => dependent_field%DECOMPOSITION%DOMAIN(mesh_component_2)%PTR% &
2855 & topology%ELEMENTS%ELEMENTS(element_number)%BASIS
2857 quadrature_scheme_2 => dependent_basis_2%QUADRATURE% &
2862 DO ns=1,dependent_basis_2%NUMBER_OF_ELEMENT_PARAMETERS
2872 pgm=quadrature_scheme_1%GAUSS_BASIS_FNS(ms,
no_part_deriv,ng)
2873 pgn=quadrature_scheme_2%GAUSS_BASIS_FNS(ns,
no_part_deriv,ng)
2876 coupling_param=equations%INTERPOLATION%MATERIALS_INTERP_POINT(field_v_variable_type)%PTR% &
2881 coupling_matrices(num_var_count)%PTR%ELEMENT_MATRIX%MATRIX(mhs,nhs) = &
2882 & coupling_matrices(num_var_count)%PTR%ELEMENT_MATRIX%MATRIX(mhs,nhs) + &
2883 & coupling_param * pgm * pgn * rwg
2900 IF(
ASSOCIATED(equations_set%ANALYTIC))
THEN 2909 DO mh=1,field_variable%NUMBER_OF_COMPONENTS
2910 mesh_component_1=field_variable%COMPONENTS(mh)%MESH_COMPONENT_NUMBER
2911 dependent_basis_1=>dependent_field%DECOMPOSITION%DOMAIN(mesh_component_1)%PTR% &
2912 & topology%ELEMENTS%ELEMENTS(element_number)%BASIS
2914 rwg=equations%INTERPOLATION%GEOMETRIC_INTERP_POINT_METRICS(field_u_variable_type)%PTR%JACOBIAN* &
2915 & quadrature_scheme_1%GAUSS_WEIGHTS(ng)
2916 DO ms=1,dependent_basis_1%NUMBER_OF_ELEMENT_PARAMETERS
2918 pgm=quadrature_scheme_1%GAUSS_BASIS_FNS(ms,
no_part_deriv,ng)
2922 x(1) = equations%INTERPOLATION%GEOMETRIC_INTERP_POINT(field_u_variable_type)%PTR%VALUES(1,1)
2923 x(2) = equations%INTERPOLATION%GEOMETRIC_INTERP_POINT(field_u_variable_type)%PTR%VALUES(2,1)
2924 IF(dependent_basis_1%NUMBER_OF_XI==3)
THEN 2925 x(3) = equations%INTERPOLATION%GEOMETRIC_INTERP_POINT(field_u_variable_type)%PTR%VALUES(3,1)
2931 fact = perm_over_vis_param / l
2934 source = -2.0_dp / l * fact * exp( arg(1) ) * exp( arg(2) )
2941 fact = 2.0_dp *
pi * perm_over_vis_param / l
2942 arg(1) = 2.0_dp *
pi * x(1) / l
2943 arg(2) = 2.0_dp *
pi * x(2) / l
2944 source = +2.0_dp * (2.0_dp *
pi / l) * fact * sin( arg(1) ) * sin( arg(2) )
2953 fact = perm_over_vis_param / l
2957 source = -3.0_dp / l * fact * exp( arg(1) ) * exp( arg(2) ) * exp( arg(3) )
2964 fact = 2.0_dp *
pi * perm_over_vis_param / l
2965 arg(1) = 2.0_dp *
pi * x(1) / l
2966 arg(2) = 2.0_dp *
pi * x(2) / l
2967 arg(3) = 2.0_dp *
pi * x(3) / l
2968 source = +3.0_dp * ( 2.0_dp *
pi / l ) * fact * sin( arg(1) ) * sin( arg(2) ) * sin( arg(3) )
2976 rhs_vector%ELEMENT_VECTOR%VECTOR(mhs)=rhs_vector%ELEMENT_VECTOR%VECTOR(mhs)+sum*rwg
2980 rhs_vector%ELEMENT_VECTOR%VECTOR(mhs)=0.0_dp
3057 IF(rhs_vector%UPDATE_VECTOR)
THEN 3064 IF( element_number == 1 )
THEN 3066 DO mh=1,field_variable%NUMBER_OF_COMPONENTS
3067 mesh_component_1 = field_variable%COMPONENTS(mh)%MESH_COMPONENT_NUMBER
3068 dependent_basis_1 => dependent_field%DECOMPOSITION%DOMAIN(mesh_component_1)%PTR% &
3069 & topology%ELEMENTS%ELEMENTS(element_number)%BASIS
3070 ndofs = ndofs + dependent_basis_1%NUMBER_OF_ELEMENT_PARAMETERS
3077 & stiffness_matrix%ELEMENT_MATRIX%MATRIX(mhs,:), &
3078 &
'("",4(X,E13.6))',
'4(4(X,E13.6))',err,error,*999)
3085 IF(dependent_field%SCALINGS%SCALING_TYPE/=field_no_scaling)
THEN 3086 CALL field_interpolationparametersscalefactorselementget(element_number,equations%INTERPOLATION% &
3087 & dependent_interp_parameters(field_var_type)%PTR,err,error,*999)
3089 DO mh=1,field_variable%NUMBER_OF_COMPONENTS
3091 mesh_component_1=field_variable%COMPONENTS(mh)%MESH_COMPONENT_NUMBER
3092 dependent_basis_1=>dependent_field%DECOMPOSITION%DOMAIN(mesh_component_1)%PTR% &
3093 & topology%ELEMENTS%ELEMENTS(element_number)%BASIS
3094 DO ms=1,dependent_basis_1%NUMBER_OF_ELEMENT_PARAMETERS
3097 IF(
ASSOCIATED(stiffness_matrix).AND.
ASSOCIATED(damping_matrix))
THEN 3098 IF(stiffness_matrix%UPDATE_MATRIX.OR.damping_matrix%UPDATE_MATRIX)
THEN 3100 DO nh=1,field_variable%NUMBER_OF_COMPONENTS
3101 mesh_component_2=field_variable%COMPONENTS(nh)%MESH_COMPONENT_NUMBER
3102 dependent_basis_2=>dependent_field%DECOMPOSITION%DOMAIN(mesh_component_2)%PTR% &
3103 & topology%ELEMENTS%ELEMENTS(element_number)%BASIS
3104 DO ns=1,dependent_basis_2%NUMBER_OF_ELEMENT_PARAMETERS
3106 IF(stiffness_matrix%UPDATE_MATRIX)
THEN 3107 stiffness_matrix%ELEMENT_MATRIX%MATRIX(mhs,nhs)=stiffness_matrix%ELEMENT_MATRIX%MATRIX(mhs,nhs)* &
3108 & equations%INTERPOLATION%DEPENDENT_INTERP_PARAMETERS(field_var_type)%PTR%SCALE_FACTORS(ms,mh)* &
3109 & equations%INTERPOLATION%DEPENDENT_INTERP_PARAMETERS(field_var_type)%PTR%SCALE_FACTORS(ns,nh)
3111 IF(damping_matrix%UPDATE_MATRIX)
THEN 3112 damping_matrix%ELEMENT_MATRIX%MATRIX(mhs,nhs)=damping_matrix%ELEMENT_MATRIX%MATRIX(mhs,nhs)* &
3113 & equations%INTERPOLATION%DEPENDENT_INTERP_PARAMETERS(field_var_type)%PTR%SCALE_FACTORS(ms,mh)* &
3114 & equations%INTERPOLATION%DEPENDENT_INTERP_PARAMETERS(field_var_type)%PTR%SCALE_FACTORS(ns,nh)
3120 IF(
ASSOCIATED(rhs_vector))
THEN 3121 IF(rhs_vector%UPDATE_VECTOR) rhs_vector%ELEMENT_VECTOR%VECTOR(mhs)=rhs_vector%ELEMENT_VECTOR%VECTOR(mhs)* &
3122 & equations%INTERPOLATION%DEPENDENT_INTERP_PARAMETERS(field_var_type)%PTR%SCALE_FACTORS(ms,mh)
3124 IF(
ASSOCIATED(source_vector))
THEN 3125 IF(source_vector%UPDATE_VECTOR) source_vector%ELEMENT_VECTOR%VECTOR(mhs)= &
3126 & source_vector%ELEMENT_VECTOR%VECTOR(mhs)* &
3127 & equations%INTERPOLATION%DEPENDENT_INTERP_PARAMETERS(field_var_type)%PTR%SCALE_FACTORS(ms,mh)
3136 local_error=
"Equations set subtype "//
trim(
number_to_vstring(equations_set_subtype,
"*",err,error))// &
3137 &
" is not valid for a Darcy equation type of a fluid mechanics equations set class." 3138 CALL flagerror(local_error,err,error,*999)
3142 CALL flagerror(
"Equations set equations is not associated.",err,error,*999)
3145 CALL flagerror(
"Equations set is not associated.",err,error,*999)
3148 exits(
"DARCY_EQUATION_FINITE_ELEMENT_CALCULATE")
3150 999 errorsexits(
"DARCY_EQUATION_FINITE_ELEMENT_CALCULATE",err,error)
3164 INTEGER(INTG),
INTENT(IN) :: elementNumber
3166 INTEGER(INTG),
INTENT(OUT) :: err
3183 INTEGER(INTG) :: faceIdx, faceNumber
3184 INTEGER(INTG) :: componentIdx, gaussIdx
3185 INTEGER(INTG) :: elementBaseDofIdx, faceNodeIdx, elementNodeIdx
3186 INTEGER(INTG) :: faceNodeDerivativeIdx, meshComponentNumber, nodeDerivativeIdx, parameterIdx
3187 INTEGER(INTG) :: faceParameterIdx, elementDofIdx, normalComponentIdx
3188 REAL(DP) :: gaussWeight, normalProjection, pressureGauss
3190 enters(
"Darcy_FiniteElementFaceIntegrate",err,error,*999)
3192 NULLIFY(decomposition)
3193 NULLIFY(decompelement)
3194 NULLIFY(dependentbasis)
3196 NULLIFY(equationsmatrices)
3199 NULLIFY(facequadraturescheme)
3200 NULLIFY(dependentinterpolatedpoint)
3201 NULLIFY(dependentinterpolationparameters)
3202 NULLIFY(geometricinterpolatedpoint)
3203 NULLIFY(geometricinterpolationparameters)
3206 equations=>equationsset%EQUATIONS
3207 IF(
ASSOCIATED(equationsset))
THEN 3208 equations=>equationsset%EQUATIONS
3209 IF(
ASSOCIATED(equations))
THEN 3210 equationsmatrices=>equations%EQUATIONS_MATRICES
3211 IF(
ASSOCIATED(equationsmatrices))
THEN 3212 rhsvector=>equationsmatrices%RHS_VECTOR
3215 CALL flagerror(
"Equations set equations is not associated.",err,error,*999)
3218 CALL flagerror(
"Equations set is not associated.",err,error,*999)
3221 IF(.NOT.
ALLOCATED(equationsset%specification))
THEN 3222 CALL flagerror(
"Equations set specification is not allocated.",err,error,*999)
3223 ELSE IF(
SIZE(equationsset%specification,1)/=3)
THEN 3224 CALL flagerror(
"Equations set specification must have three entries for a Darcy type equations set.", &
3227 SELECT CASE(equationsset%specification(3))
3234 decomposition=>dependentvariable%FIELD%DECOMPOSITION
3237 meshcomponentnumber=dependentvariable%COMPONENTS(1)%MESH_COMPONENT_NUMBER
3238 dependentbasis=>decomposition%DOMAIN(meshcomponentnumber)%PTR%TOPOLOGY%ELEMENTS% &
3239 & elements(elementnumber)%BASIS
3240 decompelement=>decomposition%TOPOLOGY%ELEMENTS%ELEMENTS(elementnumber)
3243 IF(decomposition%CALCULATE_FACES)
THEN 3245 dependentinterpolationparameters=>equations%INTERPOLATION%DEPENDENT_INTERP_PARAMETERS( &
3246 & dependentvariable%VARIABLE_TYPE)%PTR
3247 dependentinterpolatedpoint=>equations%INTERPOLATION%DEPENDENT_INTERP_POINT( &
3248 & dependentvariable%VARIABLE_TYPE)%PTR
3250 DO faceidx=1,dependentbasis%NUMBER_OF_LOCAL_FACES
3252 IF(
ALLOCATED(decompelement%ELEMENT_FACES))
THEN 3253 facenumber=decompelement%ELEMENT_FACES(faceidx)
3255 CALL flagerror(
"Decomposition element faces is not allocated.",err,error,*999)
3257 face=>decomposition%TOPOLOGY%FACES%FACES(facenumber)
3260 IF(.NOT.(face%BOUNDARY_FACE)) cycle
3261 CALL field_interpolation_parameters_face_get(field_values_set_type,facenumber, &
3262 & dependentinterpolationparameters,err,error,*999)
3263 normalcomponentidx=abs(face%XI_DIRECTION)
3264 facebasis=>decomposition%DOMAIN(meshcomponentnumber)%PTR%TOPOLOGY%FACES%FACES(facenumber)%BASIS
3267 DO gaussidx=1,facequadraturescheme%NUMBER_OF_GAUSS
3268 gaussweight=facequadraturescheme%GAUSS_WEIGHTS(gaussidx)
3271 & dependentinterpolatedpoint,err,error,*999)
3272 pressuregauss=dependentinterpolatedpoint%values(4,1)
3275 geometricinterpolationparameters=>equations%INTERPOLATION%GEOMETRIC_INTERP_PARAMETERS( &
3276 & field_u_variable_type)%PTR
3277 CALL field_interpolation_parameters_element_get(field_values_set_type,elementnumber, &
3278 & geometricinterpolationparameters,err,error,*999)
3279 geometricinterpolatedpoint=>equations%INTERPOLATION%GEOMETRIC_INTERP_POINT(field_u_variable_type)%PTR
3281 & geometricinterpolatedpoint,err,error,*999)
3283 pointmetrics=>equations%INTERPOLATION%GEOMETRIC_INTERP_POINT_METRICS(field_u_variable_type)%PTR
3286 DO componentidx=1,dependentvariable%NUMBER_OF_COMPONENTS-1
3287 normalprojection=dot_product(pointmetrics%GU(normalcomponentidx,:),pointmetrics%DX_DXI(componentidx,:))
3288 IF(face%XI_DIRECTION<0)
THEN 3289 normalprojection=-normalprojection
3293 elementbasedofidx=dependentbasis%NUMBER_OF_ELEMENT_PARAMETERS*(componentidx-1)
3294 DO facenodeidx=1,facebasis%NUMBER_OF_NODES
3295 elementnodeidx=dependentbasis%NODE_NUMBERS_IN_LOCAL_FACE(facenodeidx,faceidx)
3296 DO facenodederivativeidx=1,facebasis%NUMBER_OF_DERIVATIVES(facenodeidx)
3297 nodederivativeidx=dependentbasis%DERIVATIVE_NUMBERS_IN_LOCAL_FACE(facenodederivativeidx,facenodeidx,faceidx)
3298 parameteridx=dependentbasis%ELEMENT_PARAMETER_INDEX(nodederivativeidx,elementnodeidx)
3299 faceparameteridx=facebasis%ELEMENT_PARAMETER_INDEX(facenodederivativeidx,facenodeidx)
3300 elementdofidx=elementbasedofidx+parameteridx
3301 rhsvector%ELEMENT_VECTOR%VECTOR(elementdofidx) = rhsvector%ELEMENT_VECTOR%VECTOR(elementdofidx) - &
3302 & gaussweight*pressuregauss*normalprojection* &
3303 & facequadraturescheme%GAUSS_BASIS_FNS(faceparameteridx,
no_part_deriv,gaussidx)* &
3304 & pointmetrics%JACOBIAN
3316 exits(
"Darcy_FiniteElementFaceIntegrate")
3318 999 errorsexits(
"Darcy_FiniteElementFaceIntegrate",err,error)
3331 INTEGER(INTG),
INTENT(IN) :: specification(:)
3332 INTEGER(INTG),
INTENT(OUT) :: err
3336 INTEGER(INTG) :: subtype
3338 enters(
"Darcy_EquationsSetSpecificationSet",err,error,*999)
3340 IF(
ASSOCIATED(equationsset))
THEN 3341 IF(
SIZE(specification,1)/=3)
THEN 3342 CALL flagerror(
"Equations set specification must have three entries for a Darcy type equations set.", &
3345 subtype=specification(3)
3346 SELECT CASE(subtype)
3360 localerror=
"The third equations set specification of "//
trim(
numbertovstring(specification(3),
"*",err,error))// &
3361 &
" is not valid for a Darcy type of a fluid mechanics equations set." 3362 CALL flagerror(localerror,err,error,*999)
3365 IF(
ALLOCATED(equationsset%specification))
THEN 3366 CALL flagerror(
"Equations set specification is already allocated.",err,error,*999)
3368 ALLOCATE(equationsset%specification(3),stat=err)
3369 IF(err/=0)
CALL flagerror(
"Could not allocate equations set specification.",err,error,*999)
3373 CALL flagerror(
"Equations set is not associated.",err,error,*999)
3376 exits(
"Darcy_EquationsSetSpecificationSet")
3378 999
errors(
"Darcy_EquationsSetSpecificationSet",err,error)
3379 exits(
"Darcy_EquationsSetSpecificationSet")
3393 INTEGER(INTG),
INTENT(IN) :: problemSpecification(:)
3394 INTEGER(INTG),
INTENT(OUT) :: err
3398 INTEGER(INTG) :: problemSubtype
3400 enters(
"Darcy_ProblemSpecificationSet",err,error,*998)
3402 IF(
ASSOCIATED(problem))
THEN 3403 IF(
SIZE(problemspecification,1)==3)
THEN 3404 problemsubtype=problemspecification(3)
3405 SELECT CASE(problemsubtype)
3414 localerror=
"The third problem subtype of "//
trim(
numbertovstring(problemsubtype,
"*",err,error))// &
3415 &
" is not valid for a Darcy type of a fluid mechanics problem." 3416 CALL flagerror(localerror,err,error,*998)
3418 IF(
ALLOCATED(problem%specification))
THEN 3419 CALL flagerror(
"Problem specification is already allocated.",err,error,*998)
3421 ALLOCATE(problem%specification(3),stat=err)
3422 IF(err/=0)
CALL flagerror(
"Could not allocate problem specification.",err,error,*999)
3426 CALL flagerror(
"Darcy problem specification must have three entries.",err,error,*998)
3429 CALL flagerror(
"Problem is not associated.",err,error,*998)
3432 exits(
"Darcy_ProblemSpecificationSet")
3434 999
IF(
ALLOCATED(problem%specification))
DEALLOCATE(problem%specification)
3435 998
errors(
"Darcy_ProblemSpecificationSet",err,error)
3436 exits(
"Darcy_ProblemSpecificationSet")
3452 INTEGER(INTG),
INTENT(OUT) :: ERR
3456 TYPE(
solver_type),
POINTER :: SOLVER, SOLVER_MAT_PROPERTIES
3462 enters(
"DARCY_EQUATION_PROBLEM_SETUP",err,error,*999)
3464 NULLIFY(control_loop)
3467 NULLIFY(solver_mat_properties)
3468 NULLIFY(solver_equations)
3469 NULLIFY(solver_equations_mat_properties)
3470 IF(
ASSOCIATED(problem))
THEN 3471 IF(.NOT.
ALLOCATED(problem%SPECIFICATION))
THEN 3472 CALL flagerror(
"Problem specification is not allocated.",err,error,*999)
3473 ELSE IF(
SIZE(problem%SPECIFICATION,1)<3)
THEN 3474 CALL flagerror(
"Problem specification must have three entries for a Darcy equation problem.",err,error,*999)
3476 SELECT CASE(problem%SPECIFICATION(3))
3482 SELECT CASE(problem_setup%SETUP_TYPE)
3484 SELECT CASE(problem_setup%ACTION_TYPE)
3490 local_error=
"The action type of "//
trim(
number_to_vstring(problem_setup%ACTION_TYPE,
"*",err,error))// &
3492 &
" is invalid for a standard Darcy equation." 3493 CALL flagerror(local_error,err,error,*999)
3496 SELECT CASE(problem_setup%ACTION_TYPE)
3502 control_loop_root=>problem%CONTROL_LOOP
3506 local_error=
"The action type of "//
trim(
number_to_vstring(problem_setup%ACTION_TYPE,
"*",err,error))// &
3508 &
" is invalid for a standard Darcy equation." 3509 CALL flagerror(local_error,err,error,*999)
3513 control_loop_root=>problem%CONTROL_LOOP
3515 SELECT CASE(problem_setup%ACTION_TYPE)
3531 local_error=
"The action type of "//
trim(
number_to_vstring(problem_setup%ACTION_TYPE,
"*",err,error))// &
3533 &
" is invalid for a standard Darcy equation." 3534 CALL flagerror(local_error,err,error,*999)
3537 SELECT CASE(problem_setup%ACTION_TYPE)
3540 control_loop_root=>problem%CONTROL_LOOP
3552 control_loop_root=>problem%CONTROL_LOOP
3561 local_error=
"The action type of "//
trim(
number_to_vstring(problem_setup%ACTION_TYPE,
"*",err,error))// &
3563 &
" is invalid for a standard Darcy equation." 3564 CALL flagerror(local_error,err,error,*999)
3567 local_error=
"The setup type of "//
trim(
number_to_vstring(problem_setup%SETUP_TYPE,
"*",err,error))// &
3568 &
" is invalid for a standard Darcy equation." 3569 CALL flagerror(local_error,err,error,*999)
3576 SELECT CASE(problem_setup%SETUP_TYPE)
3578 SELECT CASE(problem_setup%ACTION_TYPE)
3584 local_error=
"The action type of "//
trim(
number_to_vstring(problem_setup%ACTION_TYPE,
"*",err,error))// &
3586 &
" is invalid for a quasistatic Darcy equation." 3587 CALL flagerror(local_error,err,error,*999)
3590 SELECT CASE(problem_setup%ACTION_TYPE)
3597 control_loop_root=>problem%CONTROL_LOOP
3601 local_error=
"The action type of "//
trim(
number_to_vstring(problem_setup%ACTION_TYPE,
"*",err,error))// &
3603 &
" is invalid for a quasistatic Darcy equation." 3604 CALL flagerror(local_error,err,error,*999)
3608 control_loop_root=>problem%CONTROL_LOOP
3610 SELECT CASE(problem_setup%ACTION_TYPE)
3626 local_error=
"The action type of "//
trim(
number_to_vstring(problem_setup%ACTION_TYPE,
"*",err,error))// &
3628 &
" is invalid for a quasistatic Darcy equation." 3629 CALL flagerror(local_error,err,error,*999)
3632 SELECT CASE(problem_setup%ACTION_TYPE)
3635 control_loop_root=>problem%CONTROL_LOOP
3647 control_loop_root=>problem%CONTROL_LOOP
3656 local_error=
"The action type of "//
trim(
number_to_vstring(problem_setup%ACTION_TYPE,
"*",err,error))// &
3658 &
" is invalid for a quasistatic Darcy equation." 3659 CALL flagerror(local_error,err,error,*999)
3662 local_error=
"The setup type of "//
trim(
number_to_vstring(problem_setup%SETUP_TYPE,
"*",err,error))// &
3663 &
" is invalid for a quasistatic Darcy equation." 3664 CALL flagerror(local_error,err,error,*999)
3671 SELECT CASE(problem_setup%SETUP_TYPE)
3673 SELECT CASE(problem_setup%ACTION_TYPE)
3679 local_error=
"The action type of "//
trim(
number_to_vstring(problem_setup%ACTION_TYPE,
"*",err,error))// &
3681 &
" is invalid for an ALE Darcy equation." 3682 CALL flagerror(local_error,err,error,*999)
3685 SELECT CASE(problem_setup%ACTION_TYPE)
3692 control_loop_root=>problem%CONTROL_LOOP
3696 local_error=
"The action type of "//
trim(
number_to_vstring(problem_setup%ACTION_TYPE,
"*",err,error))// &
3698 &
" is invalid for an ALE Darcy equation." 3699 CALL flagerror(local_error,err,error,*999)
3703 control_loop_root=>problem%CONTROL_LOOP
3705 SELECT CASE(problem_setup%ACTION_TYPE)
3726 local_error=
"The action type of "//
trim(
number_to_vstring(problem_setup%ACTION_TYPE,
"*",err,error))// &
3728 &
" is invalid for an ALE Darcy equation." 3729 CALL flagerror(local_error,err,error,*999)
3732 SELECT CASE(problem_setup%ACTION_TYPE)
3735 control_loop_root=>problem%CONTROL_LOOP
3753 control_loop_root=>problem%CONTROL_LOOP
3765 local_error=
"The action type of "//
trim(
number_to_vstring(problem_setup%ACTION_TYPE,
"*",err,error))// &
3767 &
" is invalid for an ALE Darcy equation." 3768 CALL flagerror(local_error,err,error,*999)
3771 local_error=
"The setup type of "//
trim(
number_to_vstring(problem_setup%SETUP_TYPE,
"*",err,error))// &
3772 &
" is invalid for an ALE Darcy equation." 3773 CALL flagerror(local_error,err,error,*999)
3780 SELECT CASE(problem_setup%SETUP_TYPE)
3782 SELECT CASE(problem_setup%ACTION_TYPE)
3788 local_error=
"The action type of "//
trim(
number_to_vstring(problem_setup%ACTION_TYPE,
"*",err,error))// &
3790 &
" is invalid for an ALE Darcy equation." 3791 CALL flagerror(local_error,err,error,*999)
3794 SELECT CASE(problem_setup%ACTION_TYPE)
3801 control_loop_root=>problem%CONTROL_LOOP
3805 local_error=
"The action type of "//
trim(
number_to_vstring(problem_setup%ACTION_TYPE,
"*",err,error))// &
3807 &
" is invalid for an ALE Darcy equation." 3808 CALL flagerror(local_error,err,error,*999)
3812 control_loop_root=>problem%CONTROL_LOOP
3814 SELECT CASE(problem_setup%ACTION_TYPE)
3840 local_error=
"The action type of "//
trim(
number_to_vstring(problem_setup%ACTION_TYPE,
"*",err,error))// &
3842 &
" is invalid for an ALE Darcy equation." 3843 CALL flagerror(local_error,err,error,*999)
3846 SELECT CASE(problem_setup%ACTION_TYPE)
3849 control_loop_root=>problem%CONTROL_LOOP
3867 control_loop_root=>problem%CONTROL_LOOP
3879 local_error=
"The action type of "//
trim(
number_to_vstring(problem_setup%ACTION_TYPE,
"*",err,error))// &
3881 &
" is invalid for an ALE Darcy equation." 3882 CALL flagerror(local_error,err,error,*999)
3885 local_error=
"The setup type of "//
trim(
number_to_vstring(problem_setup%SETUP_TYPE,
"*",err,error))// &
3886 &
" is invalid for an ALE Darcy equation." 3887 CALL flagerror(local_error,err,error,*999)
3894 SELECT CASE(problem_setup%SETUP_TYPE)
3896 SELECT CASE(problem_setup%ACTION_TYPE)
3902 local_error=
"The action type of "//
trim(
number_to_vstring(problem_setup%ACTION_TYPE,
"*",err,error))// &
3904 &
" is invalid for a transient Darcy fluid." 3905 CALL flagerror(local_error,err,error,*999)
3908 SELECT CASE(problem_setup%ACTION_TYPE)
3915 control_loop_root=>problem%CONTROL_LOOP
3919 local_error=
"The action type of "//
trim(
number_to_vstring(problem_setup%ACTION_TYPE,
"*",err,error))// &
3921 &
" is invalid for a transient Darcy fluid." 3922 CALL flagerror(local_error,err,error,*999)
3926 control_loop_root=>problem%CONTROL_LOOP
3928 SELECT CASE(problem_setup%ACTION_TYPE)
3947 local_error=
"The action type of "//
trim(
number_to_vstring(problem_setup%ACTION_TYPE,
"*",err,error))// &
3949 &
" is invalid for a transient Darcy fluid." 3950 CALL flagerror(local_error,err,error,*999)
3953 SELECT CASE(problem_setup%ACTION_TYPE)
3956 control_loop_root=>problem%CONTROL_LOOP
3969 control_loop_root=>problem%CONTROL_LOOP
3978 local_error=
"The action type of "//
trim(
number_to_vstring(problem_setup%ACTION_TYPE,
"*",err,error))// &
3980 &
" is invalid for a transient Darcy fluid." 3981 CALL flagerror(local_error,err,error,*999)
3984 local_error=
"The setup type of "//
trim(
number_to_vstring(problem_setup%SETUP_TYPE,
"*",err,error))// &
3985 &
" is invalid for a transient Darcy fluid." 3986 CALL flagerror(local_error,err,error,*999)
3993 local_error=
"The problem subtype of "//
trim(
number_to_vstring(problem%SPECIFICATION(3),
"*",err,error))// &
3994 &
" does not equal a standard, quasistatic or ALE Darcy equation subtype." 3995 CALL flagerror(local_error,err,error,*999)
3999 CALL flagerror(
"Problem is not associated.",err,error,*999)
4003 exits(
"DARCY_EQUATION_PROBLEM_SETUP")
4006 999 errorsexits(
"DARCY_EQUATION_PROBLEM_SETUP",err,error)
4021 INTEGER(INTG),
INTENT(OUT) :: ERR
4035 INTEGER(INTG) :: solver_matrix_idx
4038 enters(
"DARCY_EQUATION_PRE_SOLVE",err,error,*999)
4040 IF(
ASSOCIATED(control_loop))
THEN 4041 IF(
ASSOCIATED(solver))
THEN 4042 IF(
ASSOCIATED(control_loop%PROBLEM))
THEN 4043 IF(.NOT.
ALLOCATED(control_loop%PROBLEM%SPECIFICATION))
THEN 4044 CALL flagerror(
"Problem specification is not allocated.",err,error,*999)
4045 ELSE IF(
SIZE(control_loop%PROBLEM%SPECIFICATION,1)<3)
THEN 4046 CALL flagerror(
"Problem specification must have three entries for a Darcy equation problem.",err,error,*999)
4049 SELECT CASE(control_loop%PROBLEM%SPECIFICATION(3))
4066 solver_equations=>solver%SOLVER_EQUATIONS
4067 IF(
ASSOCIATED(solver_equations))
THEN 4068 solver_mapping=>solver_equations%SOLVER_MAPPING
4069 IF(
ASSOCIATED(solver_mapping))
THEN 4070 solver_matrices=>solver_equations%SOLVER_MATRICES
4071 IF(
ASSOCIATED(solver_matrices))
THEN 4072 DO solver_matrix_idx=1,solver_mapping%NUMBER_OF_SOLVER_MATRICES
4073 solver_matrix=>solver_matrices%MATRICES(solver_matrix_idx)%PTR
4074 IF(
ASSOCIATED(solver_matrix))
THEN 4075 solver_matrix%UPDATE_MATRIX=.true.
4077 CALL flagerror(
"Solver Matrix is not associated.",err,error,*999)
4081 CALL flagerror(
"Solver Matrices is not associated.",err,error,*999)
4084 CALL flagerror(
"Solver mapping is not associated.",err,error,*999)
4087 CALL flagerror(
"Solver equations is not associated.",err,error,*999)
4091 SELECT CASE(control_loop%PROBLEM%SPECIFICATION(3))
4107 NULLIFY(solver_ale_darcy)
4109 equations=>solver_ale_darcy%SOLVER_EQUATIONS%SOLVER_MAPPING%EQUATIONS_SET_TO_SOLVER_MAP(1)%EQUATIONS
4110 IF(
ASSOCIATED(equations))
THEN 4111 equations_set=>equations%EQUATIONS_SET
4112 IF(
ASSOCIATED(equations_set))
THEN 4113 IF(
ASSOCIATED(equations_set%ANALYTIC))
THEN 4146 NULLIFY(solver_ale_darcy)
4148 equations=>solver_ale_darcy%SOLVER_EQUATIONS%SOLVER_MAPPING%EQUATIONS_SET_TO_SOLVER_MAP(1)%EQUATIONS
4149 IF(
ASSOCIATED(equations))
THEN 4150 equations_set=>equations%EQUATIONS_SET
4151 IF(
ASSOCIATED(equations_set))
THEN 4152 IF(
ASSOCIATED(equations_set%ANALYTIC))
THEN 4178 local_error=
"Problem subtype "//
trim(
number_to_vstring(control_loop%PROBLEM%SPECIFICATION(3),
"*",err,error))// &
4179 &
" is not valid for a Darcy fluid type of a fluid mechanics problem class." 4180 CALL flagerror(local_error,err,error,*999)
4183 CALL flagerror(
"Problem is not associated.",err,error,*999)
4186 CALL flagerror(
"Solver is not associated.",err,error,*999)
4189 CALL flagerror(
"Control loop is not associated.",err,error,*999)
4192 exits(
"DARCY_EQUATION_PRE_SOLVE")
4194 999 errorsexits(
"DARCY_EQUATION_PRE_SOLVE",err,error)
4206 INTEGER(INTG),
INTENT(OUT) :: ERR
4213 enters(
"DARCY_CONTROL_TIME_LOOP_PRE_LOOP",err,error,*999)
4216 NULLIFY(solver_darcy)
4217 NULLIFY(control_loop_darcy)
4218 IF(.NOT.
ALLOCATED(control_loop%PROBLEM%SPECIFICATION))
THEN 4219 CALL flagerror(
"Problem specification is not allocated.",err,error,*999)
4220 ELSE IF(
SIZE(control_loop%PROBLEM%SPECIFICATION,1)<3)
THEN 4221 CALL flagerror(
"Problem specification must have three entries for a Darcy equation problem.",err,error,*999)
4223 SELECT CASE(control_loop%PROBLEM%SPECIFICATION(3))
4254 IF(control_loop%TIME_LOOP%ITERATION_NUMBER==0)
THEN 4267 exits(
"DARCY_CONTROL_TIME_LOOP_PRE_LOOP")
4269 999 errorsexits(
"DARCY_CONTROL_TIME_LOOP_PRE_LOOP",err,error)
4283 INTEGER(INTG),
INTENT(OUT) :: ERR
4286 TYPE(
field_type),
POINTER :: DEPENDENT_FIELD, GEOMETRIC_FIELD
4295 REAL(DP),
POINTER :: INITIAL_VALUES(:)
4297 INTEGER(INTG) :: FIELD_VAR_TYPE
4298 INTEGER(INTG) :: NDOFS_TO_PRINT,equations_set_idx
4301 enters(
"Darcy_PreSolveStoreReferenceData",err,error,*999)
4303 NULLIFY(solver_equations)
4304 NULLIFY(solver_mapping)
4305 NULLIFY(equations_set)
4307 IF(
ASSOCIATED(control_loop))
THEN 4308 IF(
ASSOCIATED(solver))
THEN 4310 IF(
ASSOCIATED(control_loop%PROBLEM))
THEN 4311 IF(.NOT.
ALLOCATED(control_loop%PROBLEM%SPECIFICATION))
THEN 4312 CALL flagerror(
"Problem specification is not allocated.",err,error,*999)
4313 ELSE IF(
SIZE(control_loop%PROBLEM%SPECIFICATION,1)<3)
THEN 4314 CALL flagerror(
"Problem specification must have three entries for a Darcy equation problem.",err,error,*999)
4316 SELECT CASE(control_loop%PROBLEM%SPECIFICATION(3))
4326 solver_equations=>solver%SOLVER_EQUATIONS
4327 IF(
ASSOCIATED(solver_equations))
THEN 4328 solver_mapping=>solver_equations%SOLVER_MAPPING
4329 IF(
ASSOCIATED(solver_mapping))
THEN 4330 DO equations_set_idx=1,solver_mapping%NUMBER_OF_EQUATIONS_SETS
4331 equations_set=>solver_mapping%EQUATIONS_SETS(equations_set_idx)%PTR
4332 IF(
ASSOCIATED(equations_set))
THEN 4333 IF(.NOT.
ALLOCATED(equations_set%SPECIFICATION))
THEN 4334 CALL flagerror(
"Equations set specification is not allocated.",err,error,*999)
4335 ELSE IF(
SIZE(equations_set%SPECIFICATION,1)/=3)
THEN 4336 CALL flagerror(
"Equations set specification must have three entries for a Darcy type equations set.", &
4339 SELECT CASE(equations_set%SPECIFICATION(3))
4350 dependent_field=>equations_set%DEPENDENT%DEPENDENT_FIELD
4351 geometric_field=>equations_set%GEOMETRY%GEOMETRIC_FIELD
4352 IF(
ASSOCIATED(dependent_field).AND.
ASSOCIATED(geometric_field))
THEN 4355 CALL field_parameter_sets_copy(geometric_field,field_u_variable_type, &
4356 & field_values_set_type,field_initial_values_set_type,alpha,err,error,*999)
4357 equations_mapping=>equations_set%EQUATIONS%EQUATIONS_MAPPING
4358 IF(
ASSOCIATED(equations_mapping))
THEN 4360 SELECT CASE(equations_set%SPECIFICATION(3))
4362 field_variable=>equations_mapping%LINEAR_MAPPING%EQUATIONS_MATRIX_TO_VAR_MAPS(1)%VARIABLE
4367 field_variable=>equations_mapping%DYNAMIC_MAPPING%EQUATIONS_MATRIX_TO_VAR_MAPS(1)%VARIABLE
4370 IF(
ASSOCIATED(field_variable))
THEN 4371 field_var_type=field_variable%VARIABLE_TYPE
4374 CALL field_parameter_sets_copy(dependent_field,field_var_type, &
4375 & field_values_set_type,field_initial_values_set_type,alpha,err,error,*999)
4378 NULLIFY(initial_values)
4379 CALL field_parameter_set_data_get(dependent_field,field_var_type, &
4380 & field_initial_values_set_type,initial_values,err,error,*999)
4381 ndofs_to_print =
SIZE(initial_values,1)
4384 &
'(" DEPENDENT_FIELD,FIELD_U_VARIABLE_TYPE,FIELD_INITIAL_VALUES_SET_TYPE = ",4(X,E13.6))', &
4385 &
'4(4(X,E13.6))',err,error,*999)
4386 CALL field_parameter_set_data_restore(dependent_field,field_var_type, &
4387 & field_initial_values_set_type,initial_values,err,error,*999)
4390 CALL flagerror(
"FIELD_VAR_TYPE is not associated.",err,error,*999)
4393 CALL flagerror(
"EQUATIONS_MAPPING is not associated.",err,error,*999)
4396 CALL flagerror(
"Dependent field and / or geometric field is / are not associated.",err,error,*999)
4399 local_error=
"Equations set subtype " &
4401 &
" is not valid for a Darcy equation fluid type of a fluid mechanics problem class." 4402 CALL flagerror(local_error,err,error,*999)
4405 CALL flagerror(
"Equations set is not associated.",err,error,*999)
4409 CALL flagerror(
"Solver mapping is not associated.",err,error,*999)
4412 CALL flagerror(
"Solver equations is not associated.",err,error,*999)
4415 local_error=
"Problem subtype "//
trim(
number_to_vstring(control_loop%PROBLEM%SPECIFICATION(3),
"*",err,error))// &
4416 &
" is not valid for a Darcy equation fluid type of a fluid mechanics problem class." 4417 CALL flagerror(local_error,err,error,*999)
4420 CALL flagerror(
"Problem is not associated.",err,error,*999)
4426 CALL flagerror(
"Solver is not associated.",err,error,*999)
4429 CALL flagerror(
"Control loop is not associated.",err,error,*999)
4433 exits(
"Darcy_PreSolveStoreReferenceData")
4435 999 errorsexits(
"Darcy_PreSolveStoreReferenceData",err,error)
4450 INTEGER(INTG),
INTENT(OUT) :: ERR
4462 enters(
"Darcy_PreSolveStorePreviousData",err,error,*999)
4464 NULLIFY(solver_equations)
4465 NULLIFY(solver_mapping)
4466 NULLIFY(equations_set)
4468 IF(
ASSOCIATED(control_loop))
THEN 4469 IF(
ASSOCIATED(solver))
THEN 4471 IF(
ASSOCIATED(control_loop%PROBLEM))
THEN 4472 IF(.NOT.
ALLOCATED(control_loop%PROBLEM%SPECIFICATION))
THEN 4473 CALL flagerror(
"Problem specification is not allocated.",err,error,*999)
4474 ELSE IF(
SIZE(control_loop%PROBLEM%SPECIFICATION,1)<3)
THEN 4475 CALL flagerror(
"Problem specification must have three entries for a Darcy equation problem.",err,error,*999)
4477 SELECT CASE(control_loop%PROBLEM%SPECIFICATION(3))
4487 solver_equations=>solver%SOLVER_EQUATIONS
4488 IF(
ASSOCIATED(solver_equations))
THEN 4489 solver_mapping=>solver_equations%SOLVER_MAPPING
4490 IF(
ASSOCIATED(solver_mapping))
THEN 4491 equations_set=>solver_mapping%EQUATIONS_SETS(1)%PTR
4492 IF(
ASSOCIATED(equations_set))
THEN 4493 IF(.NOT.
ALLOCATED(equations_set%SPECIFICATION))
THEN 4494 CALL flagerror(
"Equations set specification is not allocated.",err,error,*999)
4495 ELSE IF(
SIZE(equations_set%SPECIFICATION,1)/=3)
THEN 4496 CALL flagerror(
"Equations set specification must have three entries for a Darcy type equations set.", &
4499 SELECT CASE(equations_set%SPECIFICATION(3))
4510 geometric_field=>equations_set%GEOMETRY%GEOMETRIC_FIELD
4511 IF(
ASSOCIATED(geometric_field))
THEN 4514 CALL field_parameter_sets_copy(geometric_field,field_u_variable_type, &
4515 & field_values_set_type,field_previous_values_set_type,alpha,err,error,*999)
4517 CALL flagerror(
"Dependent field and / or geometric field is / are not associated.",err,error,*999)
4520 local_error=
"Equations set subtype " &
4522 &
" is not valid for a Darcy equation fluid type of a fluid mechanics problem class." 4523 CALL flagerror(local_error,err,error,*999)
4526 CALL flagerror(
"Equations set is not associated.",err,error,*999)
4529 CALL flagerror(
"Solver mapping is not associated.",err,error,*999)
4532 CALL flagerror(
"Solver equations is not associated.",err,error,*999)
4535 local_error=
"Problem subtype "//
trim(
number_to_vstring(control_loop%PROBLEM%SPECIFICATION(3),
"*",err,error))// &
4536 &
" is not valid for a Darcy equation fluid type of a fluid mechanics problem class." 4537 CALL flagerror(local_error,err,error,*999)
4540 CALL flagerror(
"Problem is not associated.",err,error,*999)
4546 CALL flagerror(
"Solver is not associated.",err,error,*999)
4549 CALL flagerror(
"Control loop is not associated.",err,error,*999)
4552 exits(
"Darcy_PreSolveStorePreviousData")
4554 999 errorsexits(
"Darcy_PreSolveStorePreviousData",err,error)
4569 INTEGER(INTG),
INTENT(OUT) :: ERR
4580 REAL(DP) :: CURRENT_TIME,TIME_INCREMENT,ALPHA
4581 REAL(DP),
POINTER :: MESH_DISPLACEMENT_VALUES(:)
4583 INTEGER(INTG) :: dof_number,NUMBER_OF_DOFS,NDOFS_TO_PRINT,loop_idx
4584 INTEGER(INTG) :: PROBLEM_SUBTYPE
4586 enters(
"DARCY_EQUATION_PRE_SOLVE_ALE_UPDATE_MESH",err,error,*999)
4588 NULLIFY(solver_ale_darcy)
4589 NULLIFY(solver_equations)
4590 NULLIFY(solver_mapping)
4591 NULLIFY(equations_set)
4593 IF(
ASSOCIATED(control_loop))
THEN 4594 control_time_loop=>control_loop
4595 DO loop_idx=1,control_loop%CONTROL_LOOP_LEVEL
4600 IF (
ASSOCIATED(control_loop%PARENT_LOOP))
THEN 4601 control_time_loop=>control_time_loop%PARENT_LOOP
4603 CALL flagerror(
"Could not find a time control loop.",err,error,*999)
4607 IF(
ASSOCIATED(solver))
THEN 4608 IF(
ASSOCIATED(control_loop%PROBLEM))
THEN 4609 IF(.NOT.
ALLOCATED(control_loop%PROBLEM%SPECIFICATION))
THEN 4610 CALL flagerror(
"Problem specification is not allocated.",err,error,*999)
4611 ELSE IF(
SIZE(control_loop%PROBLEM%SPECIFICATION,1)<3)
THEN 4612 CALL flagerror(
"Problem specification must have three entries for a Darcy equation problem.",err,error,*999)
4614 problem_subtype=control_loop%PROBLEM%SPECIFICATION(3)
4615 SELECT CASE(problem_subtype)
4626 solver_equations=>solver%SOLVER_EQUATIONS
4627 IF(
ASSOCIATED(solver_equations))
THEN 4628 solver_mapping=>solver_equations%SOLVER_MAPPING
4629 IF(
ASSOCIATED(solver_mapping))
THEN 4630 equations_set=>solver_mapping%EQUATIONS_SETS(1)%PTR
4631 IF(
ASSOCIATED(equations_set))
THEN 4632 IF(.NOT.
ALLOCATED(equations_set%SPECIFICATION))
THEN 4633 CALL flagerror(
"Equations set specification is not allocated.",err,error,*999)
4634 ELSE IF(
SIZE(equations_set%SPECIFICATION,1)/=3)
THEN 4635 CALL flagerror(
"Equations set specification must have three entries for a Darcy type equations set.", &
4638 SELECT CASE(equations_set%SPECIFICATION(3))
4650 geometric_field=>equations_set%GEOMETRY%GEOMETRIC_FIELD
4651 IF(
ASSOCIATED(geometric_field))
THEN 4653 NULLIFY(mesh_displacement_values)
4654 CALL field_parameter_set_data_get(geometric_field,field_u_variable_type, &
4655 & field_mesh_displacement_set_type,mesh_displacement_values,err,error,*999)
4657 ndofs_to_print =
SIZE(mesh_displacement_values,1)
4659 & mesh_displacement_values,
'(" MESH_DISPLACEMENT_VALUES = ",3(X,E13.6))',
'3(3(X,E13.6))', &
4663 number_of_dofs = geometric_field%VARIABLE_TYPE_MAP(field_u_variable_type)%PTR%NUMBER_OF_DOFS
4673 DO dof_number=1,number_of_dofs
4674 CALL field_parameter_set_add_local_dof(geometric_field, &
4675 & field_u_variable_type,field_values_set_type,dof_number, &
4676 & mesh_displacement_values(dof_number), &
4679 CALL field_parameter_set_update_start(geometric_field, &
4680 & field_u_variable_type, field_values_set_type,err,error,*999)
4681 CALL field_parameter_set_update_finish(geometric_field, &
4682 & field_u_variable_type, field_values_set_type,err,error,*999)
4686 alpha=1.0_dp/time_increment
4687 CALL field_parameter_sets_copy(geometric_field,field_u_variable_type, &
4688 & field_mesh_displacement_set_type,field_mesh_velocity_set_type,alpha,err,error,*999)
4689 CALL field_parameter_set_data_restore(geometric_field,field_u_variable_type, &
4690 & field_mesh_displacement_set_type,mesh_displacement_values,err,error,*999)
4692 CALL flagerror(
"Geometric field is not associated.",err,error,*999)
4695 local_error=
"Equations set subtype " &
4697 &
" is not valid for a Darcy equation fluid type of a fluid mechanics problem class." 4698 CALL flagerror(local_error,err,error,*999)
4701 CALL flagerror(
"Equations set is not associated.",err,error,*999)
4704 CALL flagerror(
"Solver mapping is not associated.",err,error,*999)
4707 CALL flagerror(
"Solver equations is not associated.",err,error,*999)
4714 &
" is not valid for a Darcy equation fluid type of a fluid mechanics problem class." 4715 CALL flagerror(local_error,err,error,*999)
4718 CALL flagerror(
"Problem is not associated.",err,error,*999)
4721 CALL flagerror(
"Solver is not associated.",err,error,*999)
4724 CALL flagerror(
"Control loop is not associated.",err,error,*999)
4727 exits(
"DARCY_EQUATION_PRE_SOLVE_ALE_UPDATE_MESH")
4729 999 errorsexits(
"DARCY_EQUATION_PRE_SOLVE_ALE_UPDATE_MESH",err,error)
4743 INTEGER(INTG),
INTENT(OUT) :: ERR
4753 TYPE(
field_type),
POINTER :: DEPENDENT_FIELD, GEOMETRIC_FIELD
4758 REAL(DP),
POINTER :: MESH_VELOCITY_VALUES(:)
4759 REAL(DP),
POINTER :: INITIAL_VALUES(:)
4760 REAL(DP),
POINTER :: DUMMY_VALUES1(:)
4761 REAL(DP) :: CURRENT_TIME,TIME_INCREMENT
4762 REAL(DP) :: PRESSURE
4764 INTEGER(INTG) :: FIELD_VAR_TYPE
4765 INTEGER(INTG) :: BOUNDARY_CONDITION_CHECK_VARIABLE
4766 INTEGER(INTG) :: dof_number,NUMBER_OF_DOFS,loop_idx
4767 INTEGER(INTG) :: NDOFS_TO_PRINT
4769 enters(
"Darcy_PreSolveUpdateBoundaryConditions",err,error,*999)
4771 IF(
ASSOCIATED(control_loop))
THEN 4772 control_time_loop=>control_loop
4773 DO loop_idx=1,control_loop%CONTROL_LOOP_LEVEL
4778 IF (
ASSOCIATED(control_loop%PARENT_LOOP))
THEN 4779 control_time_loop=>control_time_loop%PARENT_LOOP
4781 CALL flagerror(
"Could not find a time control loop.",err,error,*999)
4785 IF(
ASSOCIATED(solver))
THEN 4787 IF(
ASSOCIATED(control_loop%PROBLEM))
THEN 4788 IF(.NOT.
ALLOCATED(control_loop%PROBLEM%SPECIFICATION))
THEN 4789 CALL flagerror(
"Problem specification is not allocated.",err,error,*999)
4790 ELSE IF(
SIZE(control_loop%PROBLEM%SPECIFICATION,1)<3)
THEN 4791 CALL flagerror(
"Problem specification must have three entries for a Darcy equation problem.",err,error,*999)
4793 SELECT CASE(control_loop%PROBLEM%SPECIFICATION(3))
4803 solver_equations=>solver%SOLVER_EQUATIONS
4804 IF(
ASSOCIATED(solver_equations))
THEN 4805 solver_mapping=>solver_equations%SOLVER_MAPPING
4806 IF(
ASSOCIATED(solver_mapping))
THEN 4807 equations=>solver_mapping%EQUATIONS_SET_TO_SOLVER_MAP(1)%EQUATIONS
4808 IF(
ASSOCIATED(equations))
THEN 4809 equations_set=>equations%EQUATIONS_SET
4810 IF(
ASSOCIATED(equations_set))
THEN 4811 IF(.NOT.
ALLOCATED(equations_set%SPECIFICATION))
THEN 4812 CALL flagerror(
"Equations set specification is not allocated.",err,error,*999)
4813 ELSE IF(
SIZE(equations_set%SPECIFICATION,1)/=3)
THEN 4814 CALL flagerror(
"Equations set specification must have three entries for a Darcy type equations set.", &
4817 SELECT CASE(equations_set%SPECIFICATION(3))
4829 dependent_field=>equations_set%DEPENDENT%DEPENDENT_FIELD
4830 geometric_field=>equations_set%GEOMETRY%GEOMETRIC_FIELD
4831 IF(
ASSOCIATED(dependent_field).AND.
ASSOCIATED(geometric_field))
THEN 4832 boundary_conditions=>solver_equations%BOUNDARY_CONDITIONS
4833 IF(
ASSOCIATED(boundary_conditions))
THEN 4834 equations_mapping=>equations_set%EQUATIONS%EQUATIONS_MAPPING
4835 IF(
ASSOCIATED(equations_mapping))
THEN 4837 SELECT CASE(equations_set%SPECIFICATION(3))
4840 field_variable=>equations_mapping%LINEAR_MAPPING%EQUATIONS_MATRIX_TO_VAR_MAPS(1)%VARIABLE
4846 field_variable=>equations_mapping%DYNAMIC_MAPPING%EQUATIONS_MATRIX_TO_VAR_MAPS(1)%VARIABLE
4849 IF(
ASSOCIATED(field_variable))
THEN 4850 field_var_type=field_variable%VARIABLE_TYPE
4853 & boundary_conditions_variable,err,error,*999)
4854 IF(
ASSOCIATED(boundary_conditions_variable))
THEN 4855 NULLIFY(mesh_velocity_values)
4856 CALL field_parameter_set_data_get(geometric_field,field_u_variable_type, &
4857 & field_mesh_velocity_set_type,mesh_velocity_values,err,error,*999)
4858 NULLIFY(initial_values)
4859 CALL field_parameter_set_data_get(dependent_field,field_var_type, &
4860 & field_initial_values_set_type,initial_values,err,error,*999)
4862 NULLIFY( dummy_values1 )
4863 CALL field_parameter_set_data_get(dependent_field,field_var_type, &
4864 & field_values_set_type,dummy_values1,err,error,*999)
4865 ndofs_to_print =
SIZE(dummy_values1,1)
4867 & ndofs_to_print,dummy_values1, &
4868 &
'(" DEPENDENT_FIELD,FIELD_VAR_TYPE,FIELD_VALUES_SET_TYPE (before) = ",4(X,E13.6))', &
4869 &
'4(4(X,E13.6))',err,error,*999)
4871 number_of_dofs = dependent_field%VARIABLE_TYPE_MAP(field_var_type)%PTR%NUMBER_OF_DOFS
4872 DO dof_number=1,number_of_dofs
4873 boundary_condition_check_variable=boundary_conditions_variable% &
4874 & condition_types(dof_number)
4877 CALL field_parameter_set_update_local_dof(dependent_field, &
4878 & field_var_type,field_values_set_type,dof_number, &
4879 & initial_values(dof_number),err,error,*999)
4893 pressure = initial_values(dof_number) * (1.0_dp - exp(- current_time**2.0_dp / 0.25_dp))
4895 CALL field_parameter_set_update_local_dof(dependent_field, &
4896 & field_var_type,field_values_set_type,dof_number, &
4897 & pressure,err,error,*999)
4902 CALL field_parameter_set_update_start(dependent_field, &
4903 & field_var_type, field_values_set_type,err,error,*999)
4904 CALL field_parameter_set_update_finish(dependent_field, &
4905 & field_var_type, field_values_set_type,err,error,*999)
4907 ndofs_to_print =
SIZE(mesh_velocity_values,1)
4909 & ndofs_to_print,mesh_velocity_values, &
4910 &
'(" MESH_VELOCITY_VALUES = ",4(X,E13.6))',
'4(4(X,E13.6))',err,error,*999)
4913 ndofs_to_print =
SIZE(initial_values,1)
4915 & ndofs_to_print,initial_values, &
4916 &
'(" INITIAL_VALUES = ",4(X,E13.6))', &
4917 &
'4(4(X,E13.6))',err,error,*999)
4919 NULLIFY( dummy_values1 )
4920 CALL field_parameter_set_data_get(dependent_field,field_var_type, &
4921 & field_values_set_type,dummy_values1,err,error,*999)
4922 ndofs_to_print =
SIZE(dummy_values1,1)
4924 & ndofs_to_print,dummy_values1, &
4925 &
'(" DEPENDENT_FIELD,FIELD_VAR_TYPE,FIELD_VALUES_SET_TYPE (after) = ",4(X,E13.6))', &
4926 &
'4(4(X,E13.6))',err,error,*999)
4927 CALL field_parameter_set_data_restore(dependent_field,field_var_type, &
4928 & field_values_set_type,dummy_values1,err,error,*999)
4930 CALL field_parameter_set_data_restore(geometric_field,field_u_variable_type, &
4931 & field_mesh_velocity_set_type,mesh_velocity_values,err,error,*999)
4932 CALL field_parameter_set_data_restore(dependent_field,field_var_type, &
4933 & field_initial_values_set_type,initial_values,err,error,*999)
4935 CALL flagerror(
"Boundary condition variable is not associated.",err,error,*999)
4938 CALL field_parameter_set_update_start(equations_set%DEPENDENT%DEPENDENT_FIELD,field_var_type, &
4939 & field_values_set_type,err,error,*999)
4940 CALL field_parameter_set_update_finish(equations_set%DEPENDENT%DEPENDENT_FIELD,field_var_type, &
4941 & field_values_set_type,err,error,*999)
4944 CALL flagerror(
"FIELD_VAR_TYPE is not associated.",err,error,*999)
4947 CALL flagerror(
"EQUATIONS_MAPPING is not associated.",err,error,*999)
4950 CALL flagerror(
"Boundary conditions are not associated.",err,error,*999)
4953 CALL flagerror(
"Dependent field and/or geometric field is/are not associated.",err,error,*999)
4956 local_error=
"Equations set subtype " &
4958 &
" is not valid for a Darcy equation fluid type of a fluid mechanics problem class." 4959 CALL flagerror(local_error,err,error,*999)
4962 CALL flagerror(
"Equations set is not associated.",err,error,*999)
4965 CALL flagerror(
"Equations are not associated.",err,error,*999)
4968 CALL flagerror(
"Solver mapping is not associated.",err,error,*999)
4971 CALL flagerror(
"Solver equations are not associated.",err,error,*999)
4974 local_error=
"Problem subtype "//
trim(
number_to_vstring(control_loop%PROBLEM%SPECIFICATION(3),
"*",err,error))// &
4975 &
" is not valid for a Darcy equation fluid type of a fluid mechanics problem class." 4976 CALL flagerror(local_error,err,error,*999)
4979 CALL flagerror(
"Problem is not associated.",err,error,*999)
4985 CALL flagerror(
"Solver is not associated.",err,error,*999)
4988 CALL flagerror(
"Control loop is not associated.",err,error,*999)
4991 exits(
"Darcy_PreSolveUpdateBoundaryConditions")
4993 999
errors(
"Darcy_PreSolveUpdateBoundaryConditions",err,error)
4994 exits(
"Darcy_PreSolveUpdateBoundaryConditions")
5009 INTEGER(INTG),
INTENT(OUT) :: ERR
5013 TYPE(
solver_type),
POINTER :: SOLVER_MAT_PROPERTIES, SOLVER_ALE_DARCY
5014 TYPE(
field_type),
POINTER :: DEPENDENT_FIELD_MAT_PROPERTIES, MATERIALS_FIELD_ALE_DARCY
5016 TYPE(
solver_mapping_type),
POINTER :: SOLVER_MAPPING_MAT_PROPERTIES, SOLVER_MAPPING_ALE_DARCY
5017 TYPE(
equations_set_type),
POINTER :: EQUATIONS_SET_MAT_PROPERTIES, EQUATIONS_SET_ALE_DARCY
5020 REAL(DP),
POINTER :: DUMMY_VALUES2(:)
5022 INTEGER(INTG) :: NUMBER_OF_COMPONENTS_DEPENDENT_FIELD_MAT_PROPERTIES,NUMBER_OF_COMPONENTS_MATERIALS_FIELD_ALE_DARCY
5023 INTEGER(INTG) :: NDOFS_TO_PRINT, I
5026 enters(
"Darcy_PreSolveUpdateMatrixProperties",err,error,*999)
5028 IF(
ASSOCIATED(control_loop))
THEN 5030 NULLIFY(solver_mat_properties)
5031 NULLIFY(solver_ale_darcy)
5033 IF(
ASSOCIATED(solver))
THEN 5034 IF(
ASSOCIATED(control_loop%PROBLEM))
THEN 5035 IF(.NOT.
ALLOCATED(control_loop%PROBLEM%SPECIFICATION))
THEN 5036 CALL flagerror(
"Problem specification is not allocated.",err,error,*999)
5037 ELSE IF(
SIZE(control_loop%PROBLEM%SPECIFICATION,1)<3)
THEN 5038 CALL flagerror(
"Problem specification must have three entries for a Darcy equation problem.",err,error,*999)
5040 SELECT CASE(control_loop%PROBLEM%SPECIFICATION(3))