MATLAB Examples

# Code used in manuscript

## Section 5.1: Class polysym

```x = polysym('x'); y = polysym('y'); disp(x^2+cos(y)) ```
```x^2+cos(y) ```

Equivalence of different polysym definitions:

```clear x = polysym('x'); y = polysym('y'); %#ok<NASGU> whos ```
``` Name Size Bytes Class Attributes x 1x1 106 polysym y 1x1 106 polysym ```
```clear polysyms('x','y') whos ```
``` Name Size Bytes Class Attributes x 1x1 106 polysym y 1x1 106 polysym ```
```clear polysyms x y whos ```
``` Name Size Bytes Class Attributes x 1x1 106 polysym y 1x1 106 polysym ```

Definition of polysym vector:

```v = polysym('v',[1 2]); disp(v) ```
```v1 v2 ```
```M = polysym('M',2); disp(M) ```
```M1_1 M1_2 M2_1 M2_2 ```

Matrix operations have the same syntax as for numeric calculations:

```disp(M*v.') ```
```M1_1*v1+M1_2*v2 M2_1*v1+M2_2*v2 ```

Numbers can be entered directly or combined with a polysym object:

```x = polysym(1+1i); disp(x*(1-1i)) ```
```(1+1*I)*(1-1*I) ```

polysym variables can be input to a MATLAB function:

```f = @(y,z) z^3 + y^2*z + 1; polysyms y z disp(f(cos(y),z)) ```
```z^3+cos(y)^2*z+1 ```

Importance of the order of evaluation:

```polysyms x disp(9/16*x) disp(x*9/16) ```
```0.5625*x x*9/16 ```

Error results from an attempt to assign a polysym object to a component of a numeric array:

```y = polysym('y',2); x = zeros(2); try x(1) = y(1); %#ok<NASGU> catch error disp(error.message) end ```
```In an assignment A(I) = B, the number of elements in B and I must be the same. ```

Ways of initializing an array so it is compatible with polysym:

```x = polysym(0,size(y)); disp(x) ```
```0 0 0 0 ```
```x = y*0; disp(x) ```
```0 0 0 0 ```

## Section 5.2: Class BertiniLab

```% Solve the system \$x^2-1=0\$, \$y^2-4=0\$: poly_system = BertiniLab('variable_group',{'x','y'}, ... 'function_def',{'x^2-1';'y^2-4'}); poly_system = poly_system.solve; ```

Read solutions and match to variables

```sols = poly_system.match_solutions('real_finite_solutions'); ```

Convert output to double precision

```x = real(double(sols.x)); y = real(double(sols.y)); disp([x(:) y(:)]) ```
``` 1.0000 2.0000 1.0000 -2.0000 -1.0000 2.0000 -1.0000 -2.0000 ```

Solve , .

```polysyms z1 z2 f = @(x,y) [(29/16)*x^3 - 2*x*y; y - x^2]; poly_system = BertiniLab('function_def',f(z1,z2), ... 'variable_group',[z1 z2]); poly_system = poly_system.solve; disp(poly_system.solve_summary) ```
``` Bertini(TM) v1.4 (October 23, 2013) D.J. Bates, J.D. Hauenstein, A.J. Sommese, C.W. Wampler (using GMP v5.1.3, MPFR v3.1.2) NOTE: You have requested to use adaptive path tracking. Please make sure that you have setup the following tolerances appropriately: CoeffBound: 3.905043000000e+00, DegreeBound: 3.000000000000e+00 AMPSafetyDigits1: 1, AMPSafetyDigits2: 1, AMPMaxPrec: 1024 Tracking path 0 of 6 Finite Solution Summary NOTE: nonsingular vs singular is based on condition number and identical endpoints | Number of real solns | Number of non-real solns | Total ------------------------------------------------------------------------------------------ Non-singular | 0 | 0 | 0 Singular | 1 | 0 | 1 ------------------------------------------------------------------------------------------ Total | 1 | 0 | 1 Finite Multiplicity Summary Multiplicity | Number of real solns | Number of non-real solns ------------------------------------------------------------------------------------------ 3 | 1 | 0 ------------------------------------------------------------------------------------------ The following files may be of interest to you: main_data: A human-readable version of the solutions - main output file. raw_solutions: A list of the solutions with the corresponding path numbers. raw_data: Similar to the previous, but with the points in Bertini's homogeneous coordinates along with more information about the solutions. real_finite_solutions: A list of all real finite solutions. finite_solutions: A list of all finite solutions. nonsingular_solutions: A list of all nonsingular solutions. singular_solutions: A list of all singular solutions. ------------------------------------------------------------------------------------------ Paths Tracked: 6 Truncated infinite paths: 3 - try adjusting SecurityMaxNorm or set SecurityLevel to 1 in the input file Please see 'failed_paths' for more information about these paths. ```

Solve again with security level reset:

```poly_system.config = struct('SecurityLevel',1); poly_system = poly_system.solve; disp(poly_system.solve_summary) ```
``` Bertini(TM) v1.4 (October 23, 2013) D.J. Bates, J.D. Hauenstein, A.J. Sommese, C.W. Wampler (using GMP v5.1.3, MPFR v3.1.2) NOTE: You have requested to use adaptive path tracking. Please make sure that you have setup the following tolerances appropriately: CoeffBound: 6.725360000000e+00, DegreeBound: 3.000000000000e+00 AMPSafetyDigits1: 1, AMPSafetyDigits2: 1, AMPMaxPrec: 1024 Tracking path 0 of 6 Finite Solution Summary NOTE: nonsingular vs singular is based on condition number and identical endpoints | Number of real solns | Number of non-real solns | Total ------------------------------------------------------------------------------------------ Non-singular | 0 | 0 | 0 Singular | 1 | 0 | 1 ------------------------------------------------------------------------------------------ Total | 1 | 0 | 1 Finite Multiplicity Summary Multiplicity | Number of real solns | Number of non-real solns ------------------------------------------------------------------------------------------ 3 | 1 | 0 Infinite Solution Summary NOTE: nonsingular vs singular is based on condition number and identical endpoints | Number of real solns | Number of non-real solns | Total ------------------------------------------------------------------------------------------ Non-singular | 0 | 0 | 0 Singular | 0 | 1 | 1 ------------------------------------------------------------------------------------------ Total | 0 | 1 | 1 Infinite Multiplicity Summary Multiplicity | Number of real solns | Number of non-real solns ------------------------------------------------------------------------------------------ 3 | 0 | 1 ------------------------------------------------------------------------------------------ The following files may be of interest to you: main_data: A human-readable version of the solutions - main output file. raw_solutions: A list of the solutions with the corresponding path numbers. raw_data: Similar to the previous, but with the points in Bertini's homogeneous coordinates along with more information about the solutions. real_finite_solutions: A list of all real finite solutions. finite_solutions: A list of all finite solutions. nonsingular_solutions: A list of all nonsingular solutions. singular_solutions: A list of all singular solutions. ------------------------------------------------------------------------------------------ Paths Tracked: 6 ```

Solve using multihomogenous homotopy (needs 5 paths instead of 6):

```polysyms z1 z2 f = @(x,y) [(29/16)*x^3 - 2*x*y; y - x^2]; poly_system = BertiniLab('function_def',f(z1,z2), ... 'variable_group',{z1,z2},'config',struct('SecurityLevel',1)); poly_system = poly_system.solve; disp(poly_system.solve_summary) ```
``` Bertini(TM) v1.4 (October 23, 2013) D.J. Bates, J.D. Hauenstein, A.J. Sommese, C.W. Wampler (using GMP v5.1.3, MPFR v3.1.2) NOTE: You have requested to use adaptive path tracking. Please make sure that you have setup the following tolerances appropriately: CoeffBound: 3.856657000000e+00, DegreeBound: 4.000000000000e+00 AMPSafetyDigits1: 1, AMPSafetyDigits2: 1, AMPMaxPrec: 1024 Tracking path 0 of 5 Finite Solution Summary NOTE: nonsingular vs singular is based on condition number and identical endpoints | Number of real solns | Number of non-real solns | Total ------------------------------------------------------------------------------------------ Non-singular | 0 | 0 | 0 Singular | 1 | 0 | 1 ------------------------------------------------------------------------------------------ Total | 1 | 0 | 1 Finite Multiplicity Summary Multiplicity | Number of real solns | Number of non-real solns ------------------------------------------------------------------------------------------ 3 | 1 | 0 Infinite Solution Summary NOTE: nonsingular vs singular is based on condition number and identical endpoints | Number of real solns | Number of non-real solns | Total ------------------------------------------------------------------------------------------ Non-singular | 0 | 0 | 0 Singular | 0 | 1 | 1 ------------------------------------------------------------------------------------------ Total | 0 | 1 | 1 Infinite Multiplicity Summary Multiplicity | Number of real solns | Number of non-real solns ------------------------------------------------------------------------------------------ 2 | 0 | 1 ------------------------------------------------------------------------------------------ The following files may be of interest to you: main_data: A human-readable version of the solutions - main output file. raw_solutions: A list of the solutions with the corresponding path numbers. raw_data: Similar to the previous, but with the points in Bertini's homogeneous coordinates along with more information about the solutions. real_finite_solutions: A list of all real finite solutions. finite_solutions: A list of all finite solutions. nonsingular_solutions: A list of all nonsingular solutions. singular_solutions: A list of all singular solutions. ------------------------------------------------------------------------------------------ Paths Tracked: 5 ```

## Section 6: Customizing BertiniLab

```% A function: type('findRealRoots') ```
```function [sols,varargout] = findRealRoots(polyFun) %findRealRoots Find the real, finite, isolated roots of a polynomial %system. % % sols = findRealRoots(polyFun) inputs a set of equations (a cell array of % strings, polysym or sym object) and returns a structure with the % solutions for each variable in a field with the name of that variable. % The variable names are automatically obtained from polyFun. % % [sols,poly_system] = ... also returns a BertiniLab object, POLY_SYSTEM, % with the definition of the problem. % % Example: % sols = findRealRoots({'x^2-1'; 'y^2-1'}); % % See also: polysolve, BertiniLab % Find variables automatically polyFun = polysym(polyFun); vars = symvar(polyFun); % Solve and retrieve finite solutions. poly_system = BertiniLab('function_def',polyFun,'variable_group',vars); poly_system = solve(poly_system); sols = poly_system.match_solutions('real_finite_solutions'); % Convert solutions to double precision sols = structfun(@double,sols,'UniformOutput',false); sols = structfun(@real,sols,'UniformOutput',false); % Warn user if there is a problem. if poly_system.has_failed_paths || poly_system.has_path_crossings results = poly_system.solve_summary; istart = strfind(results,'Paths Tracked'); warning(results(istart:end)); end % Return BertiniLab object, if desired. if nargout > 1 varargout{1} = poly_system; end ```

A call to the function:

```sols = findRealRoots({'x^2-1'; 'y^2-4'}); ```

## Section 7: Subclassing and parameter homotopy

```anisotropyFcn = @(x,K) 2*[x(1)*(x(2)^2+x(3)^2); x(2)*(x(1)^2+x(3)^2); x(3)*(x(1)^2+x(2)^2) - K*x(3)]; constrFcn = @(x) deal(sum(x.^2)-1,2*x(:)); ```
```%We solve this system for a trial value of K: Kval = 0.2; N = 3; gradFcn = @(x) anisotropyFcn(x,Kval); poly_system = solveWithConstraints(gradFcn,N,constrFcn,'config',struct('SecurityLevel',1)); poly_system = solve(poly_system); disp(poly_system.solve_summary) ```
``` Bertini(TM) v1.4 (October 23, 2013) D.J. Bates, J.D. Hauenstein, A.J. Sommese, C.W. Wampler (using GMP v5.1.3, MPFR v3.1.2) NOTE: You have requested to use adaptive path tracking. Please make sure that you have setup the following tolerances appropriately: CoeffBound: 3.406454000000e+00, DegreeBound: 4.000000000000e+00 AMPSafetyDigits1: 1, AMPSafetyDigits2: 1, AMPMaxPrec: 1024 Tracking path 0 of 54 Tracking path 20 of 54 Tracking path 40 of 54 Finite Solution Summary NOTE: nonsingular vs singular is based on condition number and identical endpoints | Number of real solns | Number of non-real solns | Total ------------------------------------------------------------------------------------------ Non-singular | 26 | 0 | 26 Singular | 0 | 0 | 0 ------------------------------------------------------------------------------------------ Total | 26 | 0 | 26 Finite Multiplicity Summary Multiplicity | Number of real solns | Number of non-real solns ------------------------------------------------------------------------------------------ 1 | 26 | 0 Infinite Solution Summary NOTE: nonsingular vs singular is based on condition number and identical endpoints | Number of real solns | Number of non-real solns | Total ------------------------------------------------------------------------------------------ Non-singular | 0 | 0 | 0 Singular | 0 | 14 | 14 ------------------------------------------------------------------------------------------ Total | 0 | 14 | 14 Infinite Multiplicity Summary Multiplicity | Number of real solns | Number of non-real solns ------------------------------------------------------------------------------------------ 2 | 0 | 14 ------------------------------------------------------------------------------------------ The following files may be of interest to you: main_data: A human-readable version of the solutions - main output file. raw_solutions: A list of the solutions with the corresponding path numbers. raw_data: Similar to the previous, but with the points in Bertini's homogeneous coordinates along with more information about the solutions. real_finite_solutions: A list of all real finite solutions. finite_solutions: A list of all finite solutions. nonsingular_solutions: A list of all nonsingular solutions. singular_solutions: A list of all singular solutions. ------------------------------------------------------------------------------------------ Paths Tracked: 54 ```
```%Now we do the ab initio solve for parameter homotopy: polysyms K gradFcn = @(x) anisotropyFcn(x,K); poly_system = solveWithConstraints(gradFcn,N,constrFcn,'config', ... struct('ParameterHomotopy',1,'SecurityLevel',1),'parameter',K); % poly_system = solve(poly_system,'paramotopy.input'); ```
```ia = []; % Sometimes this method obtains the wrong number of starting points. while (numel(ia) ~= 5) poly_system = poly_system.solve('paramotopy.input'); vars = poly_system.read_solutions('nonsingular_solutions'); v = round(double(vars)*1e10); v = [sort(v(1:2,:)); v(3:4,:)]; [~,ia] = unique(abs(v.'),'rows'); end vars_reduced = vars(:,ia(:)); ```

An example of a run using parameter homotopy follows:

```poly_system.starting_points = vars_reduced; poly_system.config.ParameterHomotopy = 2; poly_system.final_parameters = 0.2; poly_system = solve(poly_system,'paramotopy.input'); disp(poly_system.solve_summary) ```
``` Bertini(TM) v1.4 (October 23, 2013) D.J. Bates, J.D. Hauenstein, A.J. Sommese, C.W. Wampler (using GMP v5.1.3, MPFR v3.1.2) NOTE: Please make sure that you have also setup the proper start points file for your parameter homotopy with the starting parameters in 'start_parameters' and the target parameters in 'final_parameters'. NOTE: You have requested to use adaptive path tracking. Please make sure that you have setup the following tolerances appropriately: CoeffBound: 1.901255000000e+00, DegreeBound: 4.000000000000e+00 AMPSafetyDigits1: 1, AMPSafetyDigits2: 1, AMPMaxPrec: 1024 Tracking path 0 of 5 Finite Solution Summary NOTE: nonsingular vs singular is based on condition number and identical endpoints | Number of real solns | Number of non-real solns | Total ------------------------------------------------------------------------------------------ Non-singular | 5 | 0 | 5 Singular | 0 | 0 | 0 ------------------------------------------------------------------------------------------ Total | 5 | 0 | 5 Finite Multiplicity Summary Multiplicity | Number of real solns | Number of non-real solns ------------------------------------------------------------------------------------------ 1 | 5 | 0 ------------------------------------------------------------------------------------------ The following files may be of interest to you: main_data: A human-readable version of the solutions - main output file. raw_solutions: A list of the solutions with the corresponding path numbers. raw_data: Similar to the previous, but with the points in Bertini's homogeneous coordinates along with more information about the solutions. real_finite_solutions: A list of all real finite solutions. finite_solutions: A list of all finite solutions. nonsingular_solutions: A list of all nonsingular solutions. singular_solutions: A list of all singular solutions. ------------------------------------------------------------------------------------------ Paths Tracked: 5 ```

Create the solution curves.

```N = 200; Kval = linspace(-1,1,N); mag = zeros(3,5,N); lam = zeros(1,5,N); problemCases = struct.empty; nCases = 0; for ii=1:N poly_system.final_parameters = Kval(ii); poly_system = solve(poly_system,'paramotopy.input'); sols = poly_system.read_solutions('raw_solutions'); if poly_system.has_failed_paths || poly_system.has_path_crossings nCases = nCases+1; problemCases(nCases).number = ii; problemCases(nCases).diagnosis = poly_system.solve_summary; end mag(:,:,ii) = sols(1:3,:); lam(:,:,ii) = sols(4,:); end if ~isempty(problemCases) warning('There were some failed paths or path crossings. See problemCases for more information.') end ```

Plot magnetization

```map = pmkmp(N,'IsoL'); %map = gray(N); for ii=1:N I = all(abs(imag(mag(:,:,ii)))<1e-12); mx = real(mag(1,I,ii)); my = real(mag(2,I,ii)); mz = real(mag(3,I,ii)); m1 = [mx -mx mx -mx mx -mx mx -mx my -my my -my my -my my -my]; m2 = [my my -my -my my my -my -my mx mx -mx -mx mx mx -mx -mx]; m3 = [mz mz mz mz -mz -mz -mz -mz mz mz mz mz -mz -mz -mz -mz]; line(m1,m2,m3,'Marker','.','LineStyle','none','Color',map(ii,:)) end colorbar caxis([min(Kval) max(Kval)]) colormap(map) view(148,26) axis off equal % Add axes and an equator line([0 0],[0 0],[0 1],'Color','k') line([0 0],[0 1],[0 0],'Color','k') line([0 1],[0 0],[0 0],'Color','k') phi = linspace(0,2*pi,200); x = cos(phi); y = sin(phi); z = zeros(size(phi)); line(x,y,z,'Color','k') % Make the colorbar a bit smaller. h_bar = findobj(gcf,'Tag','Colorbar'); initpos = get(h_bar,'Position'); set(h_bar, ... 'Position',[initpos(1)+initpos(3)*0.15 initpos(2)+initpos(4)*0.15 ... initpos(3)*0.7 initpos(4)*0.7]) ```
```%In real space there are 16 distinct curves. How many are there in complex %space? We can find out by running a positive-dimensional solve for the %system. As written, |solveWithConstraints| does not handle this, so we %will call BertiniLab directly: m = sym('m',[3 1]); syms K lambda E = K*(1-m(3)^2) + m(1)^2*m(2)^2+m(2)^2*m(3)^2+m(3)^2*m(1)^2 ... + lambda*(m(1)^2+m(2)^2+m(3)^2-1); f = [diff(E,m(1)); diff(E,m(2)); diff(E,m(3)); diff(E,lambda)]; poly_system = BertiniLab('variable_group',[m; lambda; K], ... 'function_def',f,'config',struct('TrackType',1)); poly_system = solve(poly_system); disp(poly_system.solve_summary) ```
``` Bertini(TM) v1.4 (October 23, 2013) D.J. Bates, J.D. Hauenstein, A.J. Sommese, C.W. Wampler (using GMP v5.1.3, MPFR v3.1.2) NOTE: You have requested to use adaptive path tracking. Please make sure that you have setup the following tolerances appropriately: CoeffBound: 5.535781000000e+00, DegreeBound: 3.000000000000e+00 AMPSafetyDigits1: 1, AMPSafetyDigits2: 1, AMPMaxPrec: 1024 Tracking regeneration codim 1 of 4: 3 paths to track. Tracking path 0 of 3 Sorting codimension 1 of 4: 3 paths to sort. Sorting 0 of 3 Preparing regeneration codim 2 of 4: 6 witness points to move. Moving 0 of 6 Tracking regeneration codim 2 of 4: 9 paths to track. Tracking path 0 of 9 Sorting codimension 2 of 4: 9 paths to sort. Sorting 0 of 9 Preparing regeneration codim 3 of 4: 18 witness points to move. Moving 0 of 18 Tracking regeneration codim 3 of 4: 27 paths to track. Tracking path 0 of 27 Tracking path 20 of 27 Sorting codimension 3 of 4: 27 paths to sort. Sorting 0 of 27 Sorting 20 of 27 Preparing regeneration codim 4 of 4: 24 witness points to move. Moving 0 of 24 Moving 20 of 24 Tracking regeneration codim 4 of 4: 48 paths to track. Tracking path 0 of 48 Tracking path 20 of 48 Tracking path 40 of 48 Sorting codimension 4 of 4: 48 paths to sort. Sorting 0 of 48 Sorting 20 of 48 Sorting 40 of 48 ************ Regenerative Cascade Summary ************ NOTE: nonsingular vs singular is based on rank deficiency and identical endpoints |codim| paths |witness superset| nonsingular | singular |nonsolutions| inf endpoints | other bad endpoints ---------------------------------------------------------------------------------------------------------------- | 1 | 3 | 0 | 0 | 0 | 3 | 0 | 0 | 2 | 9 | 0 | 0 | 0 | 9 | 0 | 0 | 3 | 27 | 0 | 0 | 0 | 24 | 3 | 0 | 4 | 48 | 26 | 26 | 0 | 0 | 22 | 0 ---------------------------------------------------------------------------------------------------------------- |total| 87 **************************************************** *************** Witness Set Summary **************** NOTE: nonsingular vs singular is based on rank deficiency and identical endpoints |codim| witness points | nonsingular | singular ------------------------------------------------- | 4 | 26 | 26 | 0 ------------------------------------------------- **************************************************** Calculating traces for codimension 4. Calculating 0 of 26 Calculating 20 of 26 Using monodromy to decompose codimension 4. Performing monodromy loops: 16 points left to classify Performing monodromy loops: 16 points left to classify Performing monodromy loops: 16 points left to classify Performing monodromy loops: 16 points left to classify Performing monodromy loops: 16 points left to classify Performing monodromy loops: 16 points left to classify Performing monodromy loops: 16 points left to classify Using combinatorial trace test to decompose codimension 4. ************* Witness Set Decomposition ************* | dimension | components | classified | unclassified ----------------------------------------------------- | 1 | 14 | 26 | 0 ----------------------------------------------------- ************** Decomposition by Degree ************** Dimension 1: 14 classified components ----------------------------------------------------- degree 1: 10 components degree 4: 4 components ***************************************************** ```
```BertiniClean % Delete all the files created by Bertini ```