Octopus
nl_operator.F90
Go to the documentation of this file.
1!! Copyright (C) 2002-2006 M. Marques, A. Castro, A. Rubio, G. Bertsch
2!!
3!! This program is free software; you can redistribute it and/or modify
4!! it under the terms of the GNU General Public License as published by
5!! the Free Software Foundation; either version 2, or (at your option)
6!! any later version.
7!!
8!! This program is distributed in the hope that it will be useful,
9!! but WITHOUT ANY WARRANTY; without even the implied warranty of
10!! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
11!! GNU General Public License for more details.
12!!
13!! You should have received a copy of the GNU General Public License
14!! along with this program; if not, write to the Free Software
15!! Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA
16!! 02110-1301, USA.
17!!
18
19#include "global.h"
20
24 use accel_oct_m
25 use batch_oct_m
27 use debug_oct_m
28 use global_oct_m
29 use index_oct_m
30 use iso_c_binding
31 use math_oct_m
32 use mesh_oct_m
34 use mpi_oct_m
39 use parser_oct_m
41 use space_oct_m
43 use types_oct_m
45
46 implicit none
47
48 private
49 public :: &
69
75 private
76 integer :: nri = 0
77 integer, allocatable :: imin(:)
78 integer, allocatable :: imax(:)
79 integer, allocatable :: ri(:, :)
80 end type nl_operator_index_t
81
83 type nl_operator_t
84 private
85 type(stencil_t), public :: stencil
86 type(mesh_t), pointer :: mesh => null()
87 integer, allocatable :: nn(:)
88 integer, public :: np = 0
89 ! When running in parallel mode, the next three arrays are unique on each node.
90 real(real64), allocatable, public :: w(:,:)
91
92 logical, public :: const_w = .true.
93
94 type(accel_mem_t), public :: buff_weights
95 type(accel_mem_t), public :: buff_half_weights
96
97 character(len=40) :: label
98
100 integer, public :: nri = 0
101 integer, allocatable, public :: ri(:,:)
102 integer, allocatable, public :: rimap(:)
103 integer, allocatable, public :: rimap_inv(:)
104
105 integer :: ninner = 0
106 integer :: nouter = 0
107
108 type(nl_operator_index_t) :: inner
109 type(nl_operator_index_t) :: outer
110
111 type(accel_kernel_t) :: kernel
112 type(accel_mem_t) :: buff_imin
113 type(accel_mem_t) :: buff_imax
114 type(accel_mem_t) :: buff_ri
115 type(accel_mem_t) :: buff_map
116 type(accel_mem_t) :: buff_all
117 type(accel_mem_t) :: buff_inner
118 type(accel_mem_t) :: buff_outer
119 type(accel_mem_t) :: buff_stencil
120 type(accel_mem_t) :: buff_ip_to_xyz
121 type(accel_mem_t) :: buff_xyz_to_ip
122
123 ! For multigrid solvers
124 type(nl_operator_t), public, pointer :: coarser => null()
125
126 end type nl_operator_t
127
128 integer, parameter :: &
129 OP_FORTRAN = 0, &
130 op_vec = 1, &
131 op_min = op_fortran, &
132 op_max = op_vec
133
134 integer, parameter :: &
135 OP_INVMAP = 1, &
136 op_map = 2, &
137 op_nomap = 3
138
139 integer, public, parameter :: OP_ALL = 3, op_inner = 1, op_outer = 2
140
141 interface
142 integer function op_is_available(opid, type)
143 implicit none
144 integer, intent(in) :: opid, type
145 end function op_is_available
146 end interface
147
148 integer :: dfunction_global = -1
149 integer :: zfunction_global = -1
150 integer :: function_opencl
151
152contains
153
154 ! ---------------------------------------------------------
157 subroutine nl_operator_global_init(namespace)
158 type(namespace_t), intent(in) :: namespace
159
160 integer :: default
161
163
164 !%Variable OperateDouble
165 !%Type integer
166 !%Section Execution::Optimization
167 !%Default optimized
168 !%Description
169 !% This variable selects the subroutine used to apply non-local
170 !% operators over the grid for real functions.
171 !%Option fortran 0
172 !% The standard Fortran function.
173 !%Option optimized 1
174 !% This version is optimized using vector primitives (if available).
175 !%End
176
177 !%Variable OperateComplex
178 !%Type integer
179 !%Section Execution::Optimization
180 !%Default optimized
181 !%Description
182 !% This variable selects the subroutine used to apply non-local
183 !% operators over the grid for complex functions.
184 !%Option fortran 0
185 !% The standard Fortran function.
186 !%Option optimized 1
187 !% This version is optimized using vector primitives (if available).
188 !%End
190 default = op_vec
191
192 call parse_variable(namespace, 'OperateDouble', default, dfunction_global)
193 if (.not. varinfo_valid_option('OperateDouble', dfunction_global)) call messages_input_error(namespace, 'OperateDouble')
194
195 call parse_variable(namespace, 'OperateComplex', default, zfunction_global)
196 if (.not. varinfo_valid_option('OperateComplex', zfunction_global)) call messages_input_error(namespace, 'OperateComplex')
199 if (accel_is_enabled()) then
201 !%Variable OperateAccel
202 !%Type integer
203 !%Default map
204 !%Section Execution::Optimization
205 !%Description
206 !% This variable selects the subroutine used to apply non-local
207 !% operators over the grid when an accelerator device is used.
208 !%Option invmap 1
209 !% The standard implementation ported to OpenCL.
210 !%Option map 2
211 !% A different version, more suitable for GPUs.
212 !%End
213 call parse_variable(namespace, 'OperateAccel', op_map, function_opencl)
215 call messages_obsolete_variable(namespace, 'OperateOpenCL', 'OperateAccel')
217 end if
218
220 end subroutine nl_operator_global_init
221
222 ! ---------------------------------------------------------
224 subroutine nl_operator_global_end()
225 push_sub(nl_operator_global_end)
226
228 end subroutine nl_operator_global_end
230 ! ---------------------------------------------------------
234 subroutine nl_operator_init(op, label)
235 type(nl_operator_t), intent(inout) :: op
236 character(len=*), intent(in) :: label
238 push_sub(nl_operator_init)
239
240 op%label = label
241
242 pop_sub(nl_operator_init)
243 end subroutine nl_operator_init
246 ! ---------------------------------------------------------
247 subroutine nl_operator_copy(opo, opi)
248 type(nl_operator_t), intent(inout) :: opo
249 type(nl_operator_t), target, intent(in) :: opi
250
251 push_sub(nl_operator_copy)
253 ! We cannot currently copy the GPU kernel for the nl_operator
254 assert(.not. accel_is_enabled())
255
256 call nl_operator_end(opo)
257 call nl_operator_init(opo, opi%label)
258
259 call stencil_copy(opi%stencil, opo%stencil)
260
261 opo%np = opi%np
262 opo%mesh => opi%mesh
263
264 safe_allocate_source_a(opo%nn, opi%nn)
265 safe_allocate_source_a(opo%w, opi%w)
266
267 opo%const_w = opi%const_w
268
269 opo%nri = opi%nri
270 assert(allocated(opi%ri))
271
272 safe_allocate_source_a(opo%ri, opi%ri)
273 safe_allocate_source_a(opo%rimap, opi%rimap)
274 safe_allocate_source_a(opo%rimap_inv, opi%rimap_inv)
275
276 if (opi%mesh%parallel_in_domains) then
277 opo%inner%nri = opi%inner%nri
278 safe_allocate_source_a(opo%inner%imin, opi%inner%imin)
279 safe_allocate_source_a(opo%inner%imax, opi%inner%imax)
280 safe_allocate_source_a(opo%inner%ri, opi%inner%ri)
281
282 opo%outer%nri = opi%outer%nri
283 safe_allocate_source_a(opo%outer%imin, opi%outer%imin)
284 safe_allocate_source_a(opo%outer%imax, opi%outer%imax)
285 safe_allocate_source_a(opo%outer%ri, opi%outer%ri)
286 end if
287
288 ! We create the corresponding GPU buffers
289 if (accel_is_enabled() .and. opo%const_w) then
290 call accel_create_buffer(opo%buff_weights, accel_mem_read_only, type_float, opo%stencil%size)
291 call accel_write_buffer(opo%buff_weights, opo%stencil%size, opo%w(:, 1))
292 call accel_create_buffer(opo%buff_half_weights, accel_mem_read_only, type_float, opo%stencil%size)
293 call accel_write_buffer(opo%buff_half_weights, opo%stencil%size, -m_half*opo%w(:, 1))
294 end if
295
296
297 pop_sub(nl_operator_copy)
298 end subroutine nl_operator_copy
299
300
301 ! ---------------------------------------------------------
302 subroutine nl_operator_build(space, mesh, op, np, const_w, regenerate)
303 class(space_t), intent(in) :: space
304 type(mesh_t), target, intent(in) :: mesh
305 type(nl_operator_t), intent(inout) :: op
306 integer, intent(in) :: np
307 logical, optional, intent(in) :: const_w
308 logical, optional, intent(in) :: regenerate
309
310 integer :: ii, jj, p1(space%dim), time, current
311 integer, allocatable :: st1(:), st2(:), st1r(:)
312 integer :: nn
313 integer :: ir, maxp, iinner, iouter
314 logical :: change, force_change
315 character(len=200) :: flags
316 integer, allocatable :: inner_points(:), outer_points(:), all_points(:)
317
318 push_sub(nl_operator_build)
320 if (mesh%parallel_in_domains .and. .not. const_w) then
321 call messages_experimental('Domain parallelization with curvilinear coordinates')
322 end if
323
324 assert(np > 0)
325
326 ! store values in structure
327 op%np = np
328 op%mesh => mesh
329 op%const_w = .false.
330 if (present(const_w)) op%const_w = const_w
331
332 ! allocate weights op%w
333 if (op%const_w) then
334 safe_allocate(op%w(1:op%stencil%size, 1))
335 message(1) = 'Debug: nl_operator_build: working with constant weights.'
336 call messages_info(1, debug_only=.true.)
337 else
338 safe_allocate(op%w(1:op%stencil%size, 1:op%np))
339 message(1) = 'Debug: nl_operator_build: working with non-constant weights.'
340 call messages_info(1, debug_only=.true.)
341 end if
343 ! set initially to zero
344 op%w = m_zero
345
346 ! Build lookup table
347 safe_allocate(st1(1:op%stencil%size))
348 safe_allocate(st1r(1:op%stencil%size))
349 safe_allocate(st2(1:op%stencil%size))
350
351 op%nri = 0
352 do time = 1, 2
353 st2 = 0
354 do ii = 1, np
355 p1 = 0
356 call mesh_local_index_to_coords(mesh, ii, p1)
357
358 do jj = 1, op%stencil%size
359 ! Get local index of p1 plus current stencil point.
360 st1(jj) = mesh_local_index_from_coords(mesh, p1 + op%stencil%points(:, jj))
361
362 assert(st1(jj) > 0)
363 end do
364
365 st1(1:op%stencil%size) = st1(1:op%stencil%size) - ii
366
367 change = any(st1 /= st2)
368
369 !the next is to detect when we move from a point that does not
370 !have boundary points as neighbours to one that has
371 force_change = any(st1 + ii > mesh%np) .and. all(st2 + ii - 1 <= mesh%np)
372
373 ! if the stencil changes
374 if (change .or. force_change) then
375 !store it
376 st2(:) = st1(:)
377
378 !first time, just count
379 if (time == 1) op%nri = op%nri + 1
380
381 !second time, store
382 if (time == 2) then
383 current = current + 1
384 op%ri(1:op%stencil%size, current) = st1(1:op%stencil%size)
385 end if
386 end if
387
388 if (time == 2) op%rimap(ii) = current
389
390 end do
391
392 !after counting, allocate
393 if (time == 1) then
394 safe_deallocate_a(op%ri)
395 safe_deallocate_a(op%rimap)
396 safe_deallocate_a(op%rimap_inv)
398 safe_allocate(op%ri(1:op%stencil%size, 1:op%nri))
399 safe_allocate(op%rimap(1:op%np))
400 safe_allocate(op%rimap_inv(1:op%nri + 1))
401 op%ri = 0
402 op%rimap = 0
403 op%rimap_inv = 0
404 current = 0
405
406 ! the sizes
407 if (mesh%use_curvilinear) then
408 safe_allocate(op%nn(1:op%nri))
409 ! for the moment all the sizes are the same
410 op%nn = op%stencil%size
411 end if
412 end if
413
414 end do
415
416 !the inverse mapping
417 op%rimap_inv(1) = 0
418 do jj = 1, op%np
419 op%rimap_inv(op%rimap(jj) + 1) = jj
420 end do
421 op%rimap_inv(op%nri + 1) = op%np
422
423 safe_deallocate_a(st1)
424 safe_deallocate_a(st1r)
425 safe_deallocate_a(st2)
426
427 do jj = 1, op%nri
428 nn = op%rimap_inv(jj + 1) - op%rimap_inv(jj)
429 end do
430
431 if (op%mesh%parallel_in_domains) then
432 !now build the arrays required to apply the nl_operator by parts
433
434 !count points
435 op%inner%nri = 0
436 op%outer%nri = 0
437 do ir = 1, op%nri
438 maxp = op%rimap_inv(ir + 1) + maxval(op%ri(1:op%stencil%size, ir))
439 if (maxp <= np) then
440 !inner point
441 op%inner%nri = op%inner%nri + 1
442 assert(op%inner%nri <= op%nri)
443 else
444 !outer point
445 op%outer%nri = op%outer%nri + 1
446 assert(op%outer%nri <= op%nri)
447 end if
448 end do
449
450 assert(op%inner%nri + op%outer%nri == op%nri)
451
452 if (optional_default(regenerate, .false.)) then
453 safe_deallocate_a(op%inner%imin)
454 safe_deallocate_a(op%inner%imax)
455 safe_deallocate_a(op%inner%ri)
456 safe_deallocate_a(op%outer%imin)
457 safe_deallocate_a(op%outer%imax)
458 safe_deallocate_a(op%outer%ri)
459 end if
460 safe_allocate(op%inner%imin(1:op%inner%nri + 1))
461 safe_allocate(op%inner%imax(1:op%inner%nri))
462 safe_allocate(op%inner%ri(1:op%stencil%size, 1:op%inner%nri))
463
464 safe_allocate(op%outer%imin(1:op%outer%nri + 1))
465 safe_allocate(op%outer%imax(1:op%outer%nri))
466 safe_allocate(op%outer%ri(1:op%stencil%size, 1:op%outer%nri))
467
468 !now populate the arrays
469 iinner = 0
470 iouter = 0
471 do ir = 1, op%nri
472 maxp = op%rimap_inv(ir + 1) + maxval(op%ri(1:op%stencil%size, ir))
473 if (maxp <= np) then
474 !inner point
475 iinner = iinner + 1
476 op%inner%imin(iinner) = op%rimap_inv(ir)
477 op%inner%imax(iinner) = op%rimap_inv(ir + 1)
478 op%inner%ri(1:op%stencil%size, iinner) = op%ri(1:op%stencil%size, ir)
479 else
480 !outer point
481 iouter = iouter + 1
482 op%outer%imin(iouter) = op%rimap_inv(ir)
483 op%outer%imax(iouter) = op%rimap_inv(ir + 1)
484 op%outer%ri(1:op%stencil%size, iouter) = op%ri(1:op%stencil%size, ir)
485 end if
486 end do
487
488 !verify that all points in the inner operator are actually inner
489 do ir = 1, op%inner%nri
490 do ii = op%inner%imin(ir) + 1, op%inner%imax(ir)
491 assert(all(ii + op%inner%ri(1:op%stencil%size, ir) <= mesh%np))
492 end do
493 end do
494
495 end if
496
497 if (accel_is_enabled() .and. op%const_w) then
498
499 write(flags, '(i5)') op%stencil%size
500 flags='-DNDIM=3 -DSTENCIL_SIZE='//trim(adjustl(flags))
501
502 if (op%mesh%parallel_in_domains) flags = '-DINDIRECT '//trim(flags)
503
504 select case (function_opencl)
505 case (op_invmap)
506 call accel_kernel_build(op%kernel, 'operate.cl', 'operate', flags)
507 case (op_map)
508 call accel_kernel_build(op%kernel, 'operate.cl', 'operate_map', flags)
509 end select
510
511 ! conversion to i8 needed to avoid integer overflow
512 call accel_create_buffer(op%buff_ri, accel_mem_read_only, type_integer, int(op%nri, int64)*op%stencil%size)
513 call accel_write_buffer(op%buff_ri, op%stencil%size, op%nri, op%ri)
514
515 select case (function_opencl)
516 case (op_invmap)
517 call accel_create_buffer(op%buff_imin, accel_mem_read_only, type_integer, op%nri)
518 call accel_write_buffer(op%buff_imin, op%nri, op%rimap_inv(1:))
519 call accel_create_buffer(op%buff_imax, accel_mem_read_only, type_integer, op%nri)
520 call accel_write_buffer(op%buff_imax, op%nri, op%rimap_inv(2:))
521
522 case (op_map)
523
525 call accel_write_buffer(op%buff_map, op%mesh%np, (op%rimap - 1)*op%stencil%size)
526
527 if (op%mesh%parallel_in_domains) then
528
529 safe_allocate(inner_points(1:op%mesh%np))
530 safe_allocate(outer_points(1:op%mesh%np))
531 safe_allocate(all_points(1:op%mesh%np))
532
533 op%ninner = 0
534 op%nouter = 0
535
536 do ii = 1, op%mesh%np
537 all_points(ii) = ii - 1
538 maxp = ii + maxval(op%ri(1:op%stencil%size, op%rimap(ii)))
539 if (maxp <= op%mesh%np) then
540 op%ninner = op%ninner + 1
541 inner_points(op%ninner) = ii - 1
542 else
543 op%nouter = op%nouter + 1
544 outer_points(op%nouter) = ii - 1
545 end if
546 end do
547
549 call accel_write_buffer(op%buff_all, op%mesh%np, all_points)
550
552 call accel_write_buffer(op%buff_inner, op%ninner, inner_points)
553
555 call accel_write_buffer(op%buff_outer, op%nouter, outer_points)
556
557 safe_deallocate_a(inner_points)
558 safe_deallocate_a(outer_points)
559 safe_deallocate_a(all_points)
560
561 end if
562 end select
563 end if
564
565 pop_sub(nl_operator_build)
566
567 end subroutine nl_operator_build
568
569 ! ---------------------------------------------------------
570 subroutine nl_operator_output_weights(this)
571 type(nl_operator_t), intent(inout) :: this
572
573 integer :: istencil, idir
574
576
577 write(message(1), '(3a)') 'Debug info: Finite difference weights for ', trim(this%label), '.'
578 write(message(2), '(a)') ' Spacing:'
579 do idir = 1, this%mesh%box%dim
580 write(message(2), '(a,f16.8)') trim(message(2)), this%mesh%spacing(idir)
581 end do
582 call messages_info(2, debug_only=.true.)
583
584 do istencil = 1, this%stencil%size
585 select case(this%mesh%box%dim)
586 case(1)
587 write(message(1), '(a,i3,1i4,f25.10)') ' ', istencil, this%stencil%points(1:1, istencil), this%w(istencil, 1)
588 case(2)
589 write(message(1), '(a,i3,2i4,f25.10)') ' ', istencil, this%stencil%points(1:2, istencil), this%w(istencil, 1)
590 case(3)
591 write(message(1), '(a,i3,3i4,f25.10)') ' ', istencil, this%stencil%points(1:3, istencil), this%w(istencil, 1)
592 end select
593 call messages_info(1, debug_only=.true.)
594 end do
595
597
598 end subroutine nl_operator_output_weights
599
600 ! ---------------------------------------------------------
601 subroutine nl_operator_end(op)
602 type(nl_operator_t), intent(inout) :: op
603
604 push_sub(nl_operator_end)
605
606 if (accel_is_enabled() .and. op%const_w) then
607
608 call accel_release_buffer(op%buff_ri)
609 select case (function_opencl)
610 case (op_invmap)
611 call accel_release_buffer(op%buff_imin)
612 call accel_release_buffer(op%buff_imax)
613
614 case (op_map)
615 call accel_release_buffer(op%buff_map)
616 if (op%mesh%parallel_in_domains) then
617 call accel_release_buffer(op%buff_all)
618 call accel_release_buffer(op%buff_inner)
619 call accel_release_buffer(op%buff_outer)
620 end if
621
622 case (op_nomap)
623 call accel_release_buffer(op%buff_map)
624 call accel_release_buffer(op%buff_stencil)
625 call accel_release_buffer(op%buff_xyz_to_ip)
626 call accel_release_buffer(op%buff_ip_to_xyz)
627 end select
628
629 call accel_release_buffer(op%buff_weights)
630 call accel_release_buffer(op%buff_half_weights)
631 end if
632
633 safe_deallocate_a(op%inner%imin)
634 safe_deallocate_a(op%inner%imax)
635 safe_deallocate_a(op%inner%ri)
636 safe_deallocate_a(op%outer%imin)
637 safe_deallocate_a(op%outer%imax)
638 safe_deallocate_a(op%outer%ri)
639
640 safe_deallocate_a(op%w)
641
642 safe_deallocate_a(op%ri)
643 safe_deallocate_a(op%rimap)
644 safe_deallocate_a(op%rimap_inv)
645 safe_deallocate_a(op%nn)
646
647 call stencil_end(op%stencil)
648
649 pop_sub(nl_operator_end)
650 end subroutine nl_operator_end
651
652
653 ! ---------------------------------------------------------
654 integer pure function nl_operator_get_index(op, is, ip) result(res)
655 type(nl_operator_t), intent(in) :: op
656 integer, intent(in) :: is
657 integer, intent(in) :: ip
658
659 res = ip + op%ri(is, op%rimap(ip))
660 end function nl_operator_get_index
661
662 ! ---------------------------------------------------------
663
665 type(nl_operator_t), intent(inout) :: op
666
668
669 ! Update the GPU weights
670 if (accel_is_enabled() .and. op%const_w) then
671 call accel_create_buffer(op%buff_weights, accel_mem_read_only, type_float, op%stencil%size)
672 call accel_create_buffer(op%buff_half_weights, accel_mem_read_only, type_float, op%stencil%size)
673 end if
674
677
678
679 ! ---------------------------------------------------------
680
681 subroutine nl_operator_update_gpu_buffers(op)
682 type(nl_operator_t), intent(inout) :: op
683
685
686 ! Update the GPU weights
687 if (accel_is_enabled() .and. op%const_w) then
688 call accel_write_buffer(op%buff_weights, op%stencil%size, op%w(:, 1))
689 call accel_write_buffer(op%buff_half_weights, op%stencil%size, -m_half*op%w(:, 1))
690 end if
691
693 end subroutine nl_operator_update_gpu_buffers
694
695 ! ---------------------------------------------------------
697 integer pure function nl_operator_np_zero_bc(op) result(np_bc)
698 type(nl_operator_t), intent(in) :: op
699
700 integer :: jj, ii
701
702 np_bc = 0
703 do jj = 1, op%nri
704 ii = op%rimap_inv(jj + 1) + maxval(op%ri(1:op%stencil%size, jj))
705 np_bc = max(np_bc, ii)
706 end do
707
708 end function nl_operator_np_zero_bc
709
710#include "undef.F90"
711#include "real.F90"
712#include "nl_operator_inc.F90"
713
714#include "undef.F90"
715#include "complex.F90"
716#include "nl_operator_inc.F90"
717
718end module nl_operator_oct_m
719
720!! Local Variables:
721!! mode: f90
722!! coding: utf-8
723!! End:
subroutine, public accel_kernel_build(this, file_name, kernel_name, flags)
Definition: accel.F90:1340
subroutine, public accel_release_buffer(this, async)
Definition: accel.F90:918
pure logical function, public accel_is_enabled()
Definition: accel.F90:418
integer, parameter, public accel_mem_read_only
Definition: accel.F90:195
integer pure function, public accel_max_workgroup_size()
Definition: accel.F90:1098
This module implements batches of mesh functions.
Definition: batch.F90:135
Module implementing boundary conditions in Octopus.
Definition: boundaries.F90:124
real(real64), parameter, public m_zero
Definition: global.F90:190
real(real64), parameter, public m_half
Definition: global.F90:196
This module implements the index, used for the mesh points.
Definition: index.F90:124
This module is intended to contain "only mathematical" functions and procedures.
Definition: math.F90:117
This module defines the meshes, which are used in Octopus.
Definition: mesh.F90:120
integer function, public mesh_local_index_from_coords(mesh, ix)
This function returns the local index of the point for a given vector of integer coordinates.
Definition: mesh.F90:935
subroutine, public mesh_local_index_to_coords(mesh, ip, ix)
Given a local point index, this function returns the set of integer coordinates of the point.
Definition: mesh.F90:947
subroutine, public messages_obsolete_variable(namespace, name, rep)
Definition: messages.F90:1029
character(len=256), dimension(max_lines), public message
to be output by fatal, warning
Definition: messages.F90:162
subroutine, public messages_input_error(namespace, var, details, row, column)
Definition: messages.F90:697
subroutine, public messages_experimental(name, namespace)
Definition: messages.F90:1069
subroutine, public messages_info(no_lines, iunit, debug_only, stress, all_nodes, namespace)
Definition: messages.F90:600
This module handles the communicators for the various parallelization strategies.
Definition: multicomm.F90:147
This module defines non-local operators.
subroutine, public dnl_operator_operate_diag(op, fo)
subroutine, public nl_operator_init(op, label)
initialize an instance of a non-local operator by setting the label
integer, parameter op_map
integer, parameter op_max
integer, parameter op_vec
subroutine, public dnl_operator_operate_batch(op, fi, fo, ghost_update, profile, points, factor, async)
subroutine, public dnl_operator_operate(op, fi, fo, ghost_update, profile, points)
subroutine, public nl_operator_update_gpu_buffers(op)
subroutine, public nl_operator_global_init(namespace)
initialize global settings for non-local operators
subroutine, public nl_operator_copy(opo, opi)
subroutine, public nl_operator_output_weights(this)
subroutine, public nl_operator_end(op)
integer, parameter, public op_inner
subroutine, public znl_operator_operate(op, fi, fo, ghost_update, profile, points)
subroutine, public nl_operator_global_end()
integer, parameter, public op_outer
subroutine, public nl_operator_build(space, mesh, op, np, const_w, regenerate)
subroutine, public znl_operator_operate_batch(op, fi, fo, ghost_update, profile, points, factor, async)
integer, parameter op_nomap
integer pure function, public nl_operator_np_zero_bc(op)
subroutine, public znl_operator_operate_diag(op, fo)
integer pure function, public nl_operator_get_index(op, is, ip)
subroutine, public nl_operator_allocate_gpu_buffers(op)
integer, parameter op_min
This module contains interfaces for routines in operate.c.
Definition: operate_f.F90:119
Some general things and nomenclature:
Definition: par_vec.F90:173
This module defines stencils used in Octopus.
Definition: stencil.F90:137
subroutine, public stencil_end(this)
Definition: stencil.F90:216
subroutine, public stencil_copy(input, output)
Definition: stencil.F90:198
type(type_t), public type_float
Definition: types.F90:135
type(type_t), public type_integer
Definition: types.F90:137
Describes mesh distribution to nodes.
Definition: mesh.F90:187
index type for non-local operators
data type for non local operators
int true(void)