agd.Eikonal.HFM_CUDA.graph_reverse

View Source

 1# Copyright 2020 Jean-Marie Mirebeau, University Paris-Sud, CNRS, University Paris-Saclay
 2# Distributed WITHOUT ANY WARRANTY. Licensed under the Apache License, Version 2.0, see http://www.apache.org/licenses/LICENSE-2.0
 3
 4import numpy as np
 5import cupy as cp
 6import os
 7from . import cupy_module_helper
 8from .cupy_module_helper import SetModuleConstant
 9
10def graph_reverse(fwd,fwd_weight,
11	invalid_index=None,nrev=None,blockDim=1024):
12	"""
13	Inverses a weighted graph. 
14	- invalid_index : special value to be ignored in fwd (defaults to Int_Max)
15	- nrev : expected max number of neighbors in reversed graph
16	- irev_t
17	Output : 
18	 - rev, rev_weight. !! Warning : likely not contiguous arrays. !!
19	"""
20	fwd = cp.ascontiguousarray(fwd)
21	fwd_weight = cp.ascontiguousarray(fwd_weight)
22
23	nfwd = len(fwd)
24	shape = fwd.shape[1:]
25	size = np.prod(shape)
26
27	int_t = fwd.dtype.type
28	weight_t = fwd_weight.dtype.type
29
30	if invalid_index is None: invalid_index = np.iinfo(int_t).max
31	if nrev is None: nrev = 2*nfwd # Default guess, will be increased if needed
32	for dtype in (np.int8,np.int16,np.int32,np.int64):
33		if np.iinfo(dtype).max>=nrev:
34			irev_t=dtype
35			break
36	else: raise ValueError("Impossible nrev")
37
38	rev = cp.full( (nrev,)+shape, invalid_index, dtype=int_t)
39	rev_weight = cp.zeros( (nrev,)+shape, dtype=weight_t)
40	irev = cp.zeros( (nfwd,)+shape, dtype=irev_t)
41
42	traits = {
43		'Int':int_t,
44		'weightT':weight_t,
45		'irevT':irev_t,
46		'invalid':(invalid_index,int_t),
47		'invalid_macro':True,
48	}
49
50	source = cupy_module_helper.traits_header(traits,integral_max=True)
51	cuda_path = os.path.join(os.path.dirname(os.path.realpath(__file__)),"cuda")
52	date_modified = cupy_module_helper.getmtime_max(cuda_path)
53
54	source += [
55	'#include "GraphReverse.h"',
56	f"// Date cuda code last modified : {date_modified}"]
57	cuoptions = ("-default-device", f"-I {cuda_path}") 
58
59	source="\n".join(source)
60	module = cupy_module_helper.GetModule(source,cuoptions)
61	SetModuleConstant(module,'size_tot',size,int_t)
62	SetModuleConstant(module,'nfwd',nfwd,int_t)
63	SetModuleConstant(module,'nrev',nrev,int_t)
64
65	cupy_kernel = module.get_function("GraphReverse")
66	gridDim = int(np.ceil(size/blockDim))
67
68	irev_mmax = 0
69	for i in range(nrev): # By construction, at least one reverse index is set each iter
70		irev_max = np.max(irev)
71		irev_mmax = max(irev_max,irev_mmax)
72
73		if irev_max==-1: return rev[:irev_mmax+1],rev_weight[:irev_mmax+1]
74		if irev_max==nrev: # Some vertices have large reverse multiplicities
75			return graph_reverse(fwd,fwd_weight,invalid_index,2*nrev,blockDim=blockDim)
76
77		cupy_kernel((gridDim,),(blockDim,),(fwd,rev,irev,fwd_weight,rev_weight))

def graph_reverse(fwd, fwd_weight, invalid_index=None, nrev=None, blockDim=1024): View Source

11def graph_reverse(fwd,fwd_weight,
12	invalid_index=None,nrev=None,blockDim=1024):
13	"""
14	Inverses a weighted graph. 
15	- invalid_index : special value to be ignored in fwd (defaults to Int_Max)
16	- nrev : expected max number of neighbors in reversed graph
17	- irev_t
18	Output : 
19	 - rev, rev_weight. !! Warning : likely not contiguous arrays. !!
20	"""
21	fwd = cp.ascontiguousarray(fwd)
22	fwd_weight = cp.ascontiguousarray(fwd_weight)
23
24	nfwd = len(fwd)
25	shape = fwd.shape[1:]
26	size = np.prod(shape)
27
28	int_t = fwd.dtype.type
29	weight_t = fwd_weight.dtype.type
30
31	if invalid_index is None: invalid_index = np.iinfo(int_t).max
32	if nrev is None: nrev = 2*nfwd # Default guess, will be increased if needed
33	for dtype in (np.int8,np.int16,np.int32,np.int64):
34		if np.iinfo(dtype).max>=nrev:
35			irev_t=dtype
36			break
37	else: raise ValueError("Impossible nrev")
38
39	rev = cp.full( (nrev,)+shape, invalid_index, dtype=int_t)
40	rev_weight = cp.zeros( (nrev,)+shape, dtype=weight_t)
41	irev = cp.zeros( (nfwd,)+shape, dtype=irev_t)
42
43	traits = {
44		'Int':int_t,
45		'weightT':weight_t,
46		'irevT':irev_t,
47		'invalid':(invalid_index,int_t),
48		'invalid_macro':True,
49	}
50
51	source = cupy_module_helper.traits_header(traits,integral_max=True)
52	cuda_path = os.path.join(os.path.dirname(os.path.realpath(__file__)),"cuda")
53	date_modified = cupy_module_helper.getmtime_max(cuda_path)
54
55	source += [
56	'#include "GraphReverse.h"',
57	f"// Date cuda code last modified : {date_modified}"]
58	cuoptions = ("-default-device", f"-I {cuda_path}") 
59
60	source="\n".join(source)
61	module = cupy_module_helper.GetModule(source,cuoptions)
62	SetModuleConstant(module,'size_tot',size,int_t)
63	SetModuleConstant(module,'nfwd',nfwd,int_t)
64	SetModuleConstant(module,'nrev',nrev,int_t)
65
66	cupy_kernel = module.get_function("GraphReverse")
67	gridDim = int(np.ceil(size/blockDim))
68
69	irev_mmax = 0
70	for i in range(nrev): # By construction, at least one reverse index is set each iter
71		irev_max = np.max(irev)
72		irev_mmax = max(irev_max,irev_mmax)
73
74		if irev_max==-1: return rev[:irev_mmax+1],rev_weight[:irev_mmax+1]
75		if irev_max==nrev: # Some vertices have large reverse multiplicities
76			return graph_reverse(fwd,fwd_weight,invalid_index,2*nrev,blockDim=blockDim)
77
78		cupy_kernel((gridDim,),(blockDim,),(fwd,rev,irev,fwd_weight,rev_weight))

Inverses a weighted graph.

invalid_index : special value to be ignored in fwd (defaults to Int_Max)
nrev : expected max number of neighbors in reversed graph
irev_t Output :
- rev, rev_weight. !! Warning : likely not contiguous arrays. !!