I originally thought you were working on humans. I spotted "non model organism", that changes things alot, albeit I guess humans are a non-model organism in a way. The question is what the organism. LD has been applied to every kingdom of life, except viruses. Arlequin isn't bad for relatively small to medium data sets (but tricky to get the format right - thats a separate question). If it's a genome sequence then I still suspect its plink. There are a number of recent LD programs BTW.
I did notice Colony was written in Fortran! Wow, the author must be a physicist at ZSL anyway...
If you have a defined number of SNPs - sounds small to medium data set - then I'd use Arlquein http://cmpg.unibe.ch/software/arlequin35/Arl35Methods.html, Arlquin does loads of pop gen BTW its very famous for AMOVA BTW (very cool calculation). HAPLOVIEW has been used for LD, these authors used it here, what surprised me was they avoided using Arlequin for the LD but used HAPLOVIEW. There's a plethora of recent LD algorithms, never used them so can' say, e.g. quickLD, stuff in R.
What I am definitely saying is the LD works on from raw alleles Thats the most important part of this response. I've attempted to describe the rationale below.
If this is human genomes, or large genome data I'd would recommend plink analysis performed directly on your data, as described here https://www.cog-genomics.org/plink/1.9/ld
Using any sort of matrix will remove the linkage disequilibrium signal, IF its a single matrix for all your data. If its multiple matrices you can perform a correlation statistic between the matrices - but its better just to do a straight linkage disequilibrium calculation, e.g. plink.
The reason why correlation matrices will fail is because LD needs to identify incongruence/disparity and whether that incongruence occurs within random expectations.
Basically you are making a phylogenetic tree - thats what the correlation matrix is about, but the tree assumes clonality. LD is about deviation from random expectations (equilibrium), so either a random distribution, e.g. Poison or a randomisation process is needed. What I'm attempting to say is "a tree" (single matrix) and LD are very different calculations. Whilst, technically you can do LD via matrices (plural), its better to do a LD calculation from the raw alleles.
From the comments, it appears a "tree like structure" (ancestors) is wanted to direct an LD calculation.
What you can do is subdivide the populations and perform LD on each subdivision. A tree (here using a correlation or kinship matrix) is a continuum so it's at what point do you call the the subdivision? Its doable but the information is not presented here.
What you can do is use "Structure" - this has been upgraded alot. This is a HWE, I think using Bayes and it designates how many populations (at LD) are present in your data. You can then perform a specific LD directly on each subdivision.
If you compare the output of "Structure",basically Bayesian HWE you should see a correlation with your kinship matrix. What you are inferring is a partial LD parameter that will define ancestry and whilst its doable it has never been calculated in this manner that I'm aware of: its random or it ain't and thats the criteria of significance. This is definitely how LD works in small genomes (everything below higher eukaryotes).