If $B_1$ and $B_2$ are the bases of two integer lattices $L_1$ and $L_2$, i.e.
$L_1=\{B_1n:n\in\mathbb Z^d\}$ and $L_2=\{B_2n:n\in\mathbb Z^d\}$,
is there an easy way to determine a basis for $L_1\cap L_2$? Answers of the form "Plug the matrices into a computer and ask for Hermite Normal Form, etc" are perfectly acceptable as this is a practical problem and the matrices of integers $B_1$ and $B_2$ are known, but I need some algorithmic way because the procedure will be repeated many times.
On
L1 and L2 are column vector latices
A is a 2x2 block matrix
0 is an all zeros matrix
$$
\begin{align}
A=
\begin{bmatrix}
L1 & L2\\
L1 & 0
\end{bmatrix}
\end{align}
$$
hermite normal form
$$
\begin{align}
B=hnf(A)
\end{align}
$$
$$ \begin{align} B= \begin{bmatrix} C & 0\\ D & E \end{bmatrix} \end{align} $$ B is a 2x2 block matrix where E is the latice intersection. you can get the number of columns in E by counting how many all zero columns there are in the block above.
The reason this works is because hnf produces zeros in the upper right block. These zero columns are made up of a linear combination of vectors from L1 and L2. since they add up to zero, the linear combination of L1 vectors equals the negative of the L2 vectors. This shows that these vectors are in both the L1 and L2 latice. You can recover these vectors by just looking at the L1 components, the block E matrix.
Unfortunately hnf tends to produce integer overflow errors in even small latices. However you don't need to put it in full hnf form. Any algorithm that will produce all zero columns in the top right block will work.
A set of lecture notes up at http://cseweb.ucsd.edu/classes/wi10/cse206a/lec2.pdf suggests computing the dual bases $D_1$ and $D_2$ of your lattices, then getting the HNF/orthogonalization of the concatenated matrix $[D_1\mid D_2]$ and then computing the dual of that. You'll want to be slightly careful with your algorithm for computing normal form to keep the intermediate coefficients from blowing up, but that shouldn't be too hard.