diff --git a/DESCRIPTION b/DESCRIPTION
index 42bbf678..b703a37c 100644
--- a/DESCRIPTION
+++ b/DESCRIPTION
@@ -1,8 +1,8 @@
 Package: useeior
 Type: Package
 Title: USEEIO R modeling software
-Version: 1.6.0
-Date: 2024-8-6
+Version: 1.6.1
+Date: 2024-11-4
 Authors@R: c(
   person("Ben","Young", email="ben.young@erg.com", role="aut"),
   person("Jorge","Vendries", email="jvendries@gmail.com", role="aut"),
diff --git a/NAMESPACE b/NAMESPACE
index 7f42189c..c504c06d 100644
--- a/NAMESPACE
+++ b/NAMESPACE
@@ -22,6 +22,8 @@ export(normalizeResultMatrixByTotalImpacts)
 export(plotMatrixCoefficient)
 export(printValidationResults)
 export(seeAvailableModels)
+export(testCalculationFunctions)
+export(testVisualizationFunctions)
 export(writeModelMatrices)
 export(writeModelforAPI)
 export(writeModeltoXLSX)
diff --git a/R/AdjustPrice.R b/R/AdjustPrice.R
index 3170f08a..8a473a84 100644
--- a/R/AdjustPrice.R
+++ b/R/AdjustPrice.R
@@ -19,8 +19,12 @@ adjustResultMatrixPrice <- function(matrix_name, currency_year, purchaser_price=
   }
   # Adjust price type of multiplier
   if (purchaser_price) {
-    logging::loginfo(paste("Adjusting", matrix_name, "matrix from producer to purchaser price..."))
-    mat <- adjustMultiplierPriceType(mat, currency_year, model)
+    if(is.null(model$Phi)) {
+      logging::logwarn("Model does not contain margins, purchaser price can not be calculated")
+    } else {
+      logging::loginfo(paste("Adjusting", matrix_name, "matrix from producer to purchaser price..."))
+      mat <- adjustMultiplierPriceType(mat, currency_year, model)
+    }
   } else {
     logging::loginfo(paste("Keeping", matrix_name, "matrix in producer price..."))
   }
diff --git a/R/ExternalImportFactors.R b/R/ExternalImportFactors.R
index 5eb01866..2a8f8b73 100644
--- a/R/ExternalImportFactors.R
+++ b/R/ExternalImportFactors.R
@@ -100,7 +100,7 @@ castImportFactors <- function(IFTable, model) {
 buildModelwithImportFactors <- function(model, configpaths = NULL) {
   # (see Palm et al. 2019)
 
-  logging::loginfo("Building Import A (A_m) accounting for ITA in Domestic FD.\n")
+  logging::loginfo("Building A_m (import requirements) accounting for international trade adjustment in domestic final demand.\n")
   # Re-derive import values in Use and final demand
   # _m denotes import-related structures
   model$UseTransactions_m <- model$UseTransactions - model$DomesticUseTransactions
@@ -126,5 +126,21 @@ buildModelwithImportFactors <- function(model, configpaths = NULL) {
   
   model$M_m <- M_m
   
+  model$M <- calculateMwithImportFactors(model)
+  
   return(model)
 }
+
+#' Derives an M matrix for total embodied flows from domestic and imported supply chains.
+#' @param model, An EEIO model object with model specs and crosswalk table loaded
+#' @return An M matrix of flows x sector
+calculateMwithImportFactors <- function(model) {
+  logging::loginfo("Calculating M matrix (total emissions and resource use per dollar) ...")
+  
+  # embodied flows from the use of imports by industries to make their commodities
+  # both directly (from A_m) and indirectly (by scaling it to total requirements using L_d)
+  M_mi <- model$M_m %*% model$A_m %*% model$L_d
+  
+  M <- model$M_d + M_mi
+  return(M)
+}
diff --git a/R/ValidateModel.R b/R/ValidateModel.R
index 4b41c142..bab1e4a0 100644
--- a/R/ValidateModel.R
+++ b/R/ValidateModel.R
@@ -494,3 +494,93 @@ validateHouseholdEmissions <- function(model) {
   result <- r$LCI_f[codes, names(flows)]
   all.equal(flows, result)
 }
+
+#' Test that model calculation functions are successful
+#' Includes tests for the following functions:
+#' adjustResultMatrixPrice, calculateFlowContributiontoImpact,
+#' calculateSectorContributiontoImpact, disaggregateTotalToDirectAndTier1,
+#' calculateSectorPurchasedbySectorSourcedImpact, aggregateResultMatrix,
+#' calculateMarginSectorImpacts
+#' 
+#' @param model, A fully built EEIO model object
+#' @export
+testCalculationFunctions <- function(model) {
+  target_year <- ifelse(model$specs$IOYear != 2019, 2019, 2020)
+  sector <- model$Commodities$Code_Loc[[10]]
+  indicator <- model$Indicators$meta$Name[[1]]
+  
+  matrix <- adjustResultMatrixPrice(matrix_name = "N",
+                                    currency_year = target_year,
+                                    purchaser_price = TRUE,
+                                    model)
+  if(!all(dim(model$N) == dim(matrix)) && !all(model$N == matrix)) {
+    print("Error in adjustResultMatrixPrice()")
+  }
+  
+  flow_contrib <- calculateFlowContributiontoImpact(model, sector, indicator)
+  if(!all.equal(sum(flow_contrib$contribution), 1)) {
+    print("Error in calculateFlowContributiontoImpact()")
+  }
+  
+  sector_contrib <- calculateSectorContributiontoImpact(model, sector, indicator)
+  if(!all.equal(sum(sector_contrib$contribution), 1)) {
+    print("Error in calculateSectorContributiontoImpact()")
+  }
+  
+  demand = model$DemandVectors$vectors[[1]]
+  result <- calculateSectorPurchasedbySectorSourcedImpact(y=demand, model, indicator)
+  if(model$specs$IODataSource != "stateior") {
+    # not working for 2R mode
+    agg_result <- aggregateResultMatrix(result, "Sector", model$crosswalk)
+  }
+  
+  result <- disaggregateTotalToDirectAndTier1(model, indicator)
+  
+  if(model$specs$IODataSource != "stateior") {
+    margins <- calculateMarginSectorImpacts(model)
+  }
+
+}
+
+#' Test that visualization functions are successful
+#' Includes tests for the following functions:
+#' barplotFloworImpactFractionbyRegion, barplotIndicatorScoresbySector, 
+#' heatmapSatelliteTableCoverage, heatmapSectorRanking, plotMatrixCoefficient
+#' 
+#' @param model, A fully built EEIO model object
+#' @export
+testVisualizationFunctions <- function(model) {
+  model_list <- list("model" = model)
+  loc <- model$specs$ModelRegionAcronyms[[1]]
+  indicator <- model$Indicators$meta$Name[[1]]
+  
+  fullcons <- calculateEEIOModel(model, perspective='DIRECT', demand="Consumption",
+                                 location = loc)
+  domcons <- calculateEEIOModel(model, perspective='DIRECT', demand="Consumption",
+                                location = loc,  use_domestic_requirements = TRUE)
+  barplotFloworImpactFractionbyRegion(domcons$LCIA_d,
+                                      fullcons$LCIA_d,
+                                      "Domestic Proportion of Impact")
+  ## ^^ This may not be working correctly for 2R models
+  
+  barplotIndicatorScoresbySector(model_list,
+                                 totals_by_sector_name = "GHG",
+                                 indicator_name = "Greenhouse Gases",
+                                 sector = FALSE, y_title = "")
+  
+  heatmapSatelliteTableCoverage(model, form = "Industry")
+  # ^^ not working for form = "Commodity"
+  
+  indicators <- model$Indicators$meta$Code[1:min(5, length(model$Indicators$meta$Code))]
+  
+  if(model$specs$IODataSource != "stateior") {
+    # not working for 2R models
+    heatmapSectorRanking(model, matrix = fullcons$LCIA_d, indicators,
+                         sector_to_remove = "", N_sector = 20)
+  }
+  
+  plotMatrixCoefficient(model_list, matrix_name = "D",
+                        coefficient_name = indicator,
+                        sector_to_remove = "", y_title = indicator,
+                        y_label = "Name")
+}
diff --git a/R/VisualizationFunctions.R b/R/VisualizationFunctions.R
index 7a2df283..cb757585 100644
--- a/R/VisualizationFunctions.R
+++ b/R/VisualizationFunctions.R
@@ -237,6 +237,9 @@ heatmapSatelliteTableCoverage <- function(model, form="Commodity") {
 #' @export
 heatmapSectorRanking <- function(model, matrix, indicators, sector_to_remove, N_sector, x_title = NULL,
                                  use_codes = TRUE) {
+  if(model$specs$IODataSource == "stateior") {
+    stop("heatmapSectorRanking not available for two-region models.")
+  }
   # Generate BEA sector color mapping
   mapping <- getBEASectorColorMapping(model)
   mapping$GroupName <- mapping$SectorName
diff --git a/tests/modelspecs/GAEEIOv1.0-GHG-19.yml b/inst/extdata/modelspecs/GAEEIOv1.0-GHG-19.yml
similarity index 100%
rename from tests/modelspecs/GAEEIOv1.0-GHG-19.yml
rename to inst/extdata/modelspecs/GAEEIOv1.0-GHG-19.yml
diff --git a/man/calculateMwithImportFactors.Rd b/man/calculateMwithImportFactors.Rd
new file mode 100644
index 00000000..d6725db0
--- /dev/null
+++ b/man/calculateMwithImportFactors.Rd
@@ -0,0 +1,17 @@
+% Generated by roxygen2: do not edit by hand
+% Please edit documentation in R/ExternalImportFactors.R
+\name{calculateMwithImportFactors}
+\alias{calculateMwithImportFactors}
+\title{Derives an M matrix for total embodied flows from domestic and imported supply chains.}
+\usage{
+calculateMwithImportFactors(model)
+}
+\arguments{
+\item{model, }{An EEIO model object with model specs and crosswalk table loaded}
+}
+\value{
+An M matrix of flows x sector
+}
+\description{
+Derives an M matrix for total embodied flows from domestic and imported supply chains.
+}
diff --git a/man/testCalculationFunctions.Rd b/man/testCalculationFunctions.Rd
new file mode 100644
index 00000000..b33dc432
--- /dev/null
+++ b/man/testCalculationFunctions.Rd
@@ -0,0 +1,24 @@
+% Generated by roxygen2: do not edit by hand
+% Please edit documentation in R/ValidateModel.R
+\name{testCalculationFunctions}
+\alias{testCalculationFunctions}
+\title{Test that model calculation functions are successful
+Includes tests for the following functions:
+adjustResultMatrixPrice, calculateFlowContributiontoImpact,
+calculateSectorContributiontoImpact, disaggregateTotalToDirectAndTier1,
+calculateSectorPurchasedbySectorSourcedImpact, aggregateResultMatrix,
+calculateMarginSectorImpacts}
+\usage{
+testCalculationFunctions(model)
+}
+\arguments{
+\item{model, }{A fully built EEIO model object}
+}
+\description{
+Test that model calculation functions are successful
+Includes tests for the following functions:
+adjustResultMatrixPrice, calculateFlowContributiontoImpact,
+calculateSectorContributiontoImpact, disaggregateTotalToDirectAndTier1,
+calculateSectorPurchasedbySectorSourcedImpact, aggregateResultMatrix,
+calculateMarginSectorImpacts
+}
diff --git a/man/testVisualizationFunctions.Rd b/man/testVisualizationFunctions.Rd
new file mode 100644
index 00000000..4de6fbe3
--- /dev/null
+++ b/man/testVisualizationFunctions.Rd
@@ -0,0 +1,20 @@
+% Generated by roxygen2: do not edit by hand
+% Please edit documentation in R/ValidateModel.R
+\name{testVisualizationFunctions}
+\alias{testVisualizationFunctions}
+\title{Test that visualization functions are successful
+Includes tests for the following functions:
+barplotFloworImpactFractionbyRegion, barplotIndicatorScoresbySector, 
+heatmapSatelliteTableCoverage, heatmapSectorRanking, plotMatrixCoefficient}
+\usage{
+testVisualizationFunctions(model)
+}
+\arguments{
+\item{model, }{A fully built EEIO model object}
+}
+\description{
+Test that visualization functions are successful
+Includes tests for the following functions:
+barplotFloworImpactFractionbyRegion, barplotIndicatorScoresbySector, 
+heatmapSatelliteTableCoverage, heatmapSectorRanking, plotMatrixCoefficient
+}
diff --git a/tests/test_model_build.R b/tests/test_model_build.R
index 33d0f19c..5c23720a 100644
--- a/tests/test_model_build.R
+++ b/tests/test_model_build.R
@@ -4,12 +4,14 @@
 # setwd("tests")
 library(useeior)
 # library(unittest, quietly = TRUE)
-if (!interactive()) options(warn=2, error = function() { sink(stderr()) ; traceback(3) ; q(status = 1) })
+if (!interactive()) options(warn=1, error = function() { sink(stderr()) ; traceback(3) ; q(status = 1) })
 
 ## USEEIOv2.0.1-411 Detail model with waste disaggregation
 m <- "USEEIOv2.0.1-411"
 model <- buildModel(m)
 printValidationResults(model)
+testCalculationFunctions(model)
+testVisualizationFunctions(model)
 
 ## USEEIOv2.0.1-411 Detail model with waste disaggregation (Economic only)
 m <- "USEEIOv2.0.1-411"
@@ -128,23 +130,24 @@ writeModeltoXLSX(model, ".")
 m <- "USEEIOv2.3-s-GHG-19"
 model <- buildModel(m)
 printValidationResults(model)
+testCalculationFunctions(model)
+testVisualizationFunctions(model)
 
 ## StateEEIOv1.0 Two-region Summary model
 m <- "GAEEIOv1.0-GHG-19"
-cfg <- paste0("modelspecs/", m, ".yml")
-model <- buildModel(m, configpaths = file.path(cfg))
+model <- buildModel(m)
 printValidationResults(model)
 writeModeltoXLSX(model, ".")
+testCalculationFunctions(model)
+testVisualizationFunctions(model)
 
 ## StateEEIOv1.0 Two-region Summary model (Economic only)
-model <- buildIOModel(m, configpaths = file.path(cfg))
+model <- buildIOModel(m)
 printValidationResults(model)
 writeModeltoXLSX(model, ".")
 
 ## StateEEIOv1.1 Two-region Summary model with Import Factors
-cfg <- c(paste0("modelspecs/", m, ".yml"),
-         "US_summary_import_factors_exio_2019_12sch.csv"
-         )
+cfg <- c("US_summary_import_factors_exio_2019_12sch.csv")
 model <- useeior:::initializeModel(m, configpaths = file.path(cfg))
 model$specs$Model <- "GAEEIOv1.1-GHG-19-IF"
 model$specs$ExternalImportFactors <- TRUE
@@ -157,6 +160,8 @@ model <- useeior:::loadandbuildIndicators(model)
 model <- useeior:::loadDemandVectors(model)
 model <- useeior:::constructEEIOMatrices(model, file.path(cfg))
 printValidationResults(model)
+testCalculationFunctions(model)
+testVisualizationFunctions(model)
 
 # ## StateEEIOv1.0 Two-region Summary model with "standard" Utility disaggregation
 # model <- useeior:::initializeModel(m, configpaths = file.path(cfg))