22#include < cereal/cereal.hpp>
33#include < cereal/archives/portable_binary.hpp>
44#include < cereal/types/map.hpp>
5- #include < zlib.h>
5+
6+ #include < unistd.h>
7+ #include < fcntl.h> // for open()
8+ #include < cerrno> // for errno
9+ #include < cstdio> // for perror()
610
711#include " matrix_cache.h"
812#include " moran_eigensystem.h"
@@ -42,17 +46,51 @@ static std::string store_location;
4246void init_cache (const std::string loc)
4347{
4448 store_location = loc;
49+ int fd;
50+ const std::string lockfile = store_location + " .lock" ;
51+ fd = open (lockfile.c_str (), O_WRONLY | O_CREAT, 0600 );
52+ struct flock oflock;
53+ oflock.l_type = F_WRLCK;/* Write lock*/
54+ oflock.l_whence = SEEK_SET;
55+ oflock.l_start = 0 ;
56+ oflock.l_len = 0 ;/* Lock whole file*/
57+ if (fd == -1 || fcntl (fd, F_SETLK, &oflock) == -1 )
58+ {
59+ DEBUG1 << " Couldn't acquire lock in init_cache()" ;
60+ return ;
61+ }
62+ DEBUG1 << " Acquired lock in init_cache()" ;
4563 std::ifstream in (store_location, std::ios::binary);
4664 if (in)
4765 {
4866 cereal::PortableBinaryInputArchive iarchive (in);
4967 iarchive (cache);
5068 }
69+ oflock.l_type = F_UNLCK;
70+ if (fcntl (fd, F_UNLCK, &oflock) == -1 )
71+ ERROR << " Couldn't release lock in store_cache()" ;
72+ close (fd);
73+ remove (lockfile.c_str ());
74+ DEBUG1 << " init_cache() successful" ;
5175}
5276
5377void store_cache ()
5478{
5579 DEBUG1 << " storing cache: " << store_location;
80+ int fd;
81+ const std::string lockfile = store_location + " .lock" ;
82+ fd = open (lockfile.c_str (), O_WRONLY | O_CREAT, 0600 );
83+ struct flock oflock;
84+ oflock.l_type = F_WRLCK; /* Write lock*/
85+ oflock.l_whence = SEEK_SET;
86+ oflock.l_start = 0 ;
87+ oflock.l_len = 0 ;/* Lock whole file*/
88+ if (fd == -1 || fcntl (fd, F_SETLK, &oflock) == -1 )
89+ {
90+ DEBUG1 << " Couldn't acquire lock in store_cache()" ;
91+ return ;
92+ }
93+ DEBUG1 << " Acquired lock in store_cache()" ;
5694 std::ofstream out (store_location, std::ios::binary);
5795 if (out)
5896 {
@@ -62,8 +100,13 @@ void store_cache()
62100 else
63101 {
64102 ERROR << " could not open cache file for storage" ;
65- return ;
66103 }
104+ oflock.l_type = F_UNLCK;
105+ if (fcntl (fd, F_UNLCK, &oflock) == -1 )
106+ ERROR << " Couldn't release lock in store_cache()" ;
107+ close (fd);
108+ remove (lockfile.c_str ());
109+ DEBUG1 << " store_cache() successful" ;
67110}
68111
69112typedef struct { MatrixXq coeffs; } below_coeff;
@@ -166,31 +209,6 @@ mpq_class pnkb_undist(int n, int m, int l3)
166209 return pnkb_undist_memo[key];
167210}
168211
169- template <typename Derived1, typename Derived2, typename Derived3>
170- void parallel_matmul (const Eigen::DenseBase<Derived1>& A,
171- const Eigen::DenseBase<Derived2>& B,
172- const Eigen::DenseBase<Derived3>& C,
173- Matrix<double > &dst)
174- {
175- // Compute A * diag(B) * C and store in dst.
176- // Eigen won't parallelize MatrixXq multiplication.
177- int m = A.rows ();
178- int n = C.cols ();
179- DEBUG1 << m << " " << n;
180- dst.resize (m, n);
181- MatrixXq tmp (m, n);
182- tmp.setZero ();
183- #pragma omp parallel for collapse(2)
184- for (int i = 0 ; i < m; ++i)
185- for (int j = 0 ; j < n; ++j)
186- for (int k = 0 ; k < B.size (); ++k)
187- tmp (i, j) += A (i, k) * B (k) * C (k, j);
188- #pragma omp parallel for collapse(2)
189- for (int i = 0 ; i < m; ++i)
190- for (int j = 0 ; j < n; ++j)
191- dst (i, j) = tmp (i, j).get_d ();
192- }
193-
194212MatrixCache& cached_matrices (const int n)
195213{
196214 if (cache.count (n) == 0 )
0 commit comments