1- """Utility functions used in CompEcon
1+ """
2+ Utility functions used in CompEcon
23
34Based routines found in the CompEcon toolbox by Miranda and Fackler.
45
89and Finance, MIT Press, 2002.
910
1011"""
11-
1212from functools import reduce
13-
1413import numpy as np
15- from numba import njit
1614
1715
1816def ckron (* arrays ):
19- """Repeatedly applies the np.kron function to an arbitrary number of
17+ """
18+ Repeatedly applies the np.kron function to an arbitrary number of
2019 input arrays
2120
2221 Parameters
@@ -43,7 +42,8 @@ def ckron(*arrays):
4342
4443
4544def gridmake (* arrays ):
46- """Expands one or more vectors (or matrices) into a matrix where rows span the
45+ """
46+ Expands one or more vectors (or matrices) into a matrix where rows span the
4747 cartesian product of combinations of the input arrays. Each column of the
4848 input arrays will correspond to one column of the output matrix.
4949
@@ -78,12 +78,13 @@ def gridmake(*arrays):
7878 out = _gridmake2 (out , arr )
7979
8080 return out
81- raise NotImplementedError ("Come back here" )
81+ else :
82+ raise NotImplementedError ("Come back here" )
8283
8384
84- @njit (cache = True , fastmath = True )
8585def _gridmake2 (x1 , x2 ):
86- """Expands two vectors (or matrices) into a matrix where rows span the
86+ """
87+ Expands two vectors (or matrices) into a matrix where rows span the
8788 cartesian product of combinations of the input arrays. Each column of the
8889 input arrays will correspond to one column of the output matrix.
8990
@@ -112,23 +113,11 @@ def _gridmake2(x1, x2):
112113
113114 """
114115 if x1 .ndim == 1 and x2 .ndim == 1 :
115- n1 = x1 .shape [0 ]
116- n2 = x2 .shape [0 ]
117- out = np .empty ((n1 * n2 , 2 ), dtype = x1 .dtype )
118- for i in range (n2 ):
119- for j in range (n1 ):
120- out [i * n1 + j , 0 ] = x1 [j ]
121- out [i * n1 + j , 1 ] = x2 [i ]
122- return out
123- if x1 .ndim > 1 and x2 .ndim == 1 :
124- n1 = x1 .shape [0 ]
125- n2 = x2 .shape [0 ]
126- ncol1 = x1 .shape [1 ]
127- out = np .empty ((n1 * n2 , ncol1 + 1 ), dtype = x1 .dtype )
128- for i in range (n2 ):
129- for j in range (n1 ):
130- for k in range (ncol1 ):
131- out [i * n1 + j , k ] = x1 [j , k ]
132- out [i * n1 + j , ncol1 ] = x2 [i ]
133- return out
134- raise NotImplementedError ("Come back here" )
116+ return np .column_stack ([np .tile (x1 , x2 .shape [0 ]),
117+ np .repeat (x2 , x1 .shape [0 ])])
118+ elif x1 .ndim > 1 and x2 .ndim == 1 :
119+ first = np .tile (x1 , (x2 .shape [0 ], 1 ))
120+ second = np .repeat (x2 , x1 .shape [0 ])
121+ return np .column_stack ([first , second ])
122+ else :
123+ raise NotImplementedError ("Come back here" )
0 commit comments