##### section 3 - basic usage ############################################################ ## ----load---------------------------------------------------------------- library("m2r") ## ----basic-call---------------------------------------------------------- m2("1 + 1") ## ----m2-not-parsed------------------------------------------------------- m2("1.2") ## ----m2-parse-m2--------------------------------------------------------- m2_parse(m2("1.2")) ## ----persistence--------------------------------------------------------- m2("a = 1") m2("a") ## ----m2-helpers---------------------------------------------------------- m2_ls() m2_exists(c("a", "b")) m2_getwd() ## ----ring-creation------------------------------------------------------- (R <- ring("t", "x", "y", "z", coefring = "QQ")) ## ----ideal-creation------------------------------------------------------ (I <- ideal("t^4 - x", "t^3 - y", "t^2 - z")) ## ----gb-computation-1---------------------------------------------------- gb(I) ## ----gb-computation-2---------------------------------------------------- gb("t^4 - x", "t^3 - y", "t^2 - z") ## ----standard-evaluation------------------------------------------------- polys <- c("t^4 - x", "t^3 - y", "t^2 - z") gb_(polys, ring = R) ## ----radical------------------------------------------------------------- QQx <- ring("x", coefring = "QQ") I <- ideal("x^2") radical(I) ## ----saturation---------------------------------------------------------- I <- ideal("(x-1) x (x+1)") J <- ideal("x") saturate(I, J) ## ----primary-decomposition-1--------------------------------------------- QQxyz <- ring("x", "y", "z", coefring = "QQ") I <- ideal("x z", "y z") (ideal_list <- primary_decomposition(I)) ## ----primary-decomposition-2--------------------------------------------- dimension(ideal_list) ## ----arithmetic---------------------------------------------------------- I <- ideal("x", "y") J <- ideal("z") I + J I * J I == J ## ----arithmetic-combo---------------------------------------------------- is.radical <- function (ideal) ideal == radical(ideal) is.radical(I) ## ----primary-decomposition-3--------------------------------------------- library("magrittr") ideal("x z", "y z") %>% primary_decomposition %>% dimension ## ----factor-n------------------------------------------------------------ (x <- 2^5 * 3^4 * 5^3 * 7^2 * 11^1) factor_n(x) gmp::factorize(x) ## ----factor-poly--------------------------------------------------------- QQxy <- ring("x", "y", coefring = "QQ") factor_poly("x^4 - y^4") ## ----snf----------------------------------------------------------------- M <- matrix(c( 2, 4, 4, -6, 6, 12, 10, -4, -16 ), nrow = 3, byrow = TRUE) mats <- snf(M); P <- mats$P; D <- mats$D; Q <- mats$Q P %*% M %*% Q solve(P) %*% D %*% solve(Q) det(P) ##### section 4 - applications ############################################################ ## ----gb-elimination-1---------------------------------------------------- use_ring(QQxyz) I <- ideal("x + y + z", "x^2 + y^2 + z^2 - 9", "x^2 + y^2 - z^2") (grobner_basis <- gb(I)) ## ----gb-elimination-2---------------------------------------------------- extract_unipoly <- function(mpolyList) Filter(is.unipoly, mpolyList)[[1]] which_unipoly <- function(mpolyList) which(sapply(mpolyList, is.unipoly)) solve_gb <- function(gb) { poly <- extract_unipoly(gb) elim_var <- vars(poly) solns <- solve_unipoly(poly, real_only = TRUE) if (length(gb) == 1) { return(structure(t(t(solns)), .Dimnames = list(NULL, elim_var))) } gb <- structure( gb[-which_unipoly(gb)[1], drop = FALSE], class = "mpolyList" ) new_systems <- lapply(solns, function(soln) plug(gb, elim_var, soln)) low_solns_list <- lapply(new_systems, solve_gb) lower_var_names <- colnames(low_solns_list[[1]]) Map(cbind, solns, low_solns_list) %>% do.call("rbind", .) %>% structure(.Dimnames = list(NULL, c(elim_var, lower_var_names))) } ## ----gb-elimination-3---------------------------------------------------- (solns <- solve_gb(grobner_basis)) %>% structure(.Dimnames = list(paste("Soln", 1:4, ":"), c("z", "y", "x"))) ## ----gb-elimination-4---------------------------------------------------- f <- as.function(grobner_basis, varorder = c("z", "y", "x"), vector = TRUE) apply(solns, 1, f) %>% apply(2, round, digits = 14) %>% structure(.Dimnames = list(paste("Eqn", 5:7, ":"), paste("Soln", 1:4))) ## ----numerical-strat----------------------------------------------------- resid <- function(v) { x <- v[1]; y <- v[2]; z <- v[3] (x + y + z)^2 + (x^2 + y^2 + z^2 - 9)^2 + (x^2 + y^2 - z^2)^2 } optim(c(x = 0, y = 0, z = 0), resid)$par ## ----indep-ideal--------------------------------------------------------- ring("p00", "p01", "p10", "p11", coefring = "QQ") indep_ideal <- ideal( "p00 - (p00 + p01) (p00 + p10)", "p01 - (p00 + p01) (p01 + p11)", "p10 - (p10 + p11) (p00 + p10)", "p11 - (p10 + p11) (p01 + p11)", "p00 + p01 + p10 + p11 - 1" ) gb(indep_ideal) ## ----dimension----------------------------------------------------------- dimension(indep_ideal) ##### section 5 - internals ############################################################ ## ----data-structures-1--------------------------------------------------- str(R) ## ----data-structures-2--------------------------------------------------- m2_name(R) m2_meta(R) %>% str ## ----data-structures-3--------------------------------------------------- P str(P) ## ----parse-function------------------------------------------------------ m2_parse_function.m2_ideal <- function(x) { m2_structure( m2_name = "", m2_class = "m2_ideal", m2_meta = list(rmap = x[[1]]) ) } ## ----pointer-1----------------------------------------------------------- gb.(I) ## ------------------------------------------------------------------------ use_ring(QQxyz) (I. <- ideal.("x + y + z", "x^2 + y^2 + z^2 - 9", "x^2 + y^2 - z^2")) gb(I.) ##### section 6 - connecting to Macaulay2 ############################################################ ## ----m2-start-cloud------------------------------------------------------ stop_m2() start_m2(cloud = TRUE) m2("1+1")