# Introduction of JijModeling

This section introduces the first time users of JijModeling to mathematical modeling with JijModeling.

## Mathematical Modeling

### Mathematical Model

A mathematical model is a mathematical expression of a problem. For example, let's look at a mathematical model for the problem of "selecting the smallest element from the elements of a one-dimensional list $\mathbf{d} = [d_0, d_1, \dots, d_{N-1}]$.

We aim to minimize equation (1) with $N$ binary variables $\mathbf{x} = (x_0, x_1, \dots, x_{N-1})$ expressed in equation (3). If we only minimize equation (1), then $\mathbf{x} = \mathbf{0}$. However, the restriction in equation (2) requires that only one of $x_0, x_1, \dots, x_{N-1}$ takes one and the rest are 0. Therefore, by minimizing equation (1) while satisfying equation (2), only the binary variable corresponding to the smallest element of $\mathbf{d}$ will be 1. From the above, we can see that this mathematical model can express the problem of "selecting the smallest element from a one-dimensional list.

### Mathematical model and Instances

In the above, we have not specified a specific value for $\mathbf{d}$. It could be [10, -3, 6.6, 0.1] or some other value. When solving the problem, we must determine the specific value of $\mathbf{d}$.

Such a specific value for solving the problem is called an `instance`

.
It is helpful to know the word instance because JijModeling separates mathematical models and instances.

In the above example, the mathematical model and instances are related as shown in the table below.

Mathematical model | $\min \sum_n d_n x_n,~\mathrm{s.t.} \sum_n x_n =1, x_n \in \{0, 1\}$ |
---|---|

Instance | $\mathbf{d} = [0.5, 1, -1, 1, 3]$ |

You can read instances as the specific data of the problem to be solved.

## Mathematical modeling by JijModeling

Let's start modeling with JijModeling.

As mentioned above, JijModeling separately represents "mathematical models" and "instances. So, first, we build a mathematical model that does not depend on instances. After that, we assign instances to a mathematical model.

First, we create the problem we want to solve as follows

`import jijmodeling as jm`

# create problem instance

problem = jm.Problem('my_first_problem')

After this, we will add objective functions and constraints to this `problem`

object.

### Variables for instances

Let's define variables for the instances needed in the optimization problem (e.g., distances between cities in a traveling salesman problem).

A variable $A$ can be defined as follows.

`import jijmodeling as jm`

# Define variable A

A = jm.Placeholder('A')

If you want to define a variable like $B_{ij}$, which is a two-dimensional array, use the following.

`# Define a 2-dimensional variable sequence B`

B = jm.Placeholder('B', dim=2)

We can declare only that $B_{ij}$ is a 2-dimensional list, without specifying the number of elements in $B_{ij}$. By defining it in this way, we can have the flexibility to change the number of elements in $B_{ij}$ on the instance side without having to modify this part of the definition.

If you want to pull information on the number of elements from this 2-dimensional array, use `shape`

as follows.

`# Define a 2-dimensional variable sequence B`

B = jm.Placeholder('B', dim=2)

# Extract the number of rows N in B

N = B.shape[0].

# Extract the number of columns M of B

M = B.shape[1]

You can specify the dimensions of the variables by `dim`

and get information on the number of elements in each dimensional variable from `shape`

.

If you want to specify the number of elements in advance, define the following

`# Define the number of rows C, columns D`

C = jm.Placeholder('C')

D = jm.Placeholder('D')

# Define 2-dimensional variable column E with rows C and columns D

E = jm.Placeholder('E', shape=(C, D))

Describe the number of elements in a tuple in the `shape`

argument of `jm.Placeholder`

.

### Decision variables

The decision variable used for combinatorial optimization is the binary variable $x = \{ 0, 1\}$ can be defined as follows.

`import jijmodeling as jm`

# Define a binary variable x

x = jm.Binary('x')

This allows us to define a single binary variable, `x`

.

By also using the `Placeholder`

from earlier, we can also define a three-dimensional binary variable sequence $y$.

`# Define variables F, G, H`

F = jm.Placeholder('F')

G = jm.Placeholder('G')

H = jm.Placeholder('H')

# Define a 3-dimensional binary variable sequence y for F x G x H

y = jm.Binary('y', shape=(F, G, H))

Just as you define instance variables like a multidimensional tensor with `jm.Placeholder`

, you can define a multidimensional decision variable tensor with `jm.Binary`

by specifying the number of elements in a tuple for `shape`

.

We can also use an integer variable. You can define an integer variable that satisfies $L \leq a \leq U$ as follows.

`# Define lower L and upper U`

L = jm.Placeholder('L')

U = jm.Placeholder('U')

# Define an integer variable a

a = jm.LogEncInteger('a', lower=L, upper=U)

Integer variables for multidimensional lists are also supported. To define an integer variable that satisfies $L_{ij} \leq a_{ij} \leq U_{ij} \leq (\forall i, j)$, do the following.

`# Define the number of elements in the integer variable`

I = jm.Placeholder('I')

J = jm.Placeholder('J')

# Define lower L and upper U

L = jm.Placeholder('L', shape=(I, J))

U = jm.Placeholder('U', shape=(I, J))

# Define an integer variable a

a = jm.LogEncInteger('a', lower=L, upper=U, shape=(I, J))

For such a definition, $L_{ij}, a_{ij}, U_{ij}$ must be multidimensional array of the same dimension and same shape.

If all integer variables in a multidimensional array have the same upper $u$ and lower $\ell$, such as $\ell \leq b_{ij} \leq u$, you can define the following.

`ell = jm.Placeholder('ell')`

u = jm.Placeholder('u')

I = jm.Placeholder('I')

J = jm.Placeholder('J')

b = jm.LogEncInteger('b', lower=ell, upper=u, shape=(I, J))

### Variables for indexing

When they are represented in the form of (multidimensional) arrays, such as the variable $d_n$ for an instance or the decision variable $x_n$ for optimization, we need an index that specifies how many elements they are. Here is how this is defined.

The index $n$ used for a one-dimensional binary variable $x_n$ or an instance variable $d_n$ is defined as follows using `jm.Element`

.

`# Define variable N`

N = jm.Placeholder('N')

# define subscript n with 0≤n<N

n = jm.Element('n', (0, N))

where the `(0, N)`

tuple represents the possible range of the index `n`

, as in $0 \leq n < N$.

Note that $N$ is not included when written as `(0, N)`

. `(0, N)`

means $\{0, 1, ..., N-1\}$.

Any set can be used for the part of `jm.Element`

that specifies the range. For example

`# Define a variable M`

M = jm.Placeholder('M', dim=1)

# define an index m satisfying m ∈ M

m = jm.Element('m', M)

and then write $M = \{ 2, 3, 5, 7, 11\}$ as an instance, you have defined $m \in \{ 2, 3, 5, 7, 11\}$.

In addition, a two-dimensional list or similar can be used for the part that specifies the range. For example

`# Define a 2-dimensional variable E`

E = jm.Placeholder('E', dim=2)

# define an index e satisfying e[0] ∈ E[0], e[1] ∈ E[1].

e = jm.Element('e', E)

and give $E = [[1, 2], [3, 4], [5, 6]]$ as instances, so that `e[0]`

is 1, 3, 5 and `e[1]`

will act as an index that takes values only for 2, 4, and 6. In this way, for example, for a set of edges $E$ in a graph, $(uv) \in E$, we can define the vertex $(uv)$ that forms that edge.

## Basic Summation operation

Here is how to implement sums of sums appearing in a mathematical model using JijModeling.

Let's start with the basic

is an implementation of an operation such as It is written as follows

`# Define variables to be used`

N = jm.Placeholder('N')

x = jm.Binary('x', shape=(N))

n = jm.Element('i', (0, N))

# Calculate the sum

sum_x = jm.Sum(n, x[n])

The desired sum operation can be implemented by entering the subscripts to be summed as the first argument of `jm.Sum`

and the formula to be summed as the second argument.

$\mathrm{jm.Sum(\overbrace{n}^{subscript}, \underbrace{x[n]}_{operand})}$

Substituting the index defined in `jm.Element`

for the index to be summed, the sum of all possible ranges for that index is calculated. In the present case, the operation is to find the sum of $n = 0, 1, \dots, N-1$.

Such a simple sum can be written in abbreviated form as follows

`# Define variables to be used`

N = jm.Placeholder('N')

x = jm.Binary('x', shape=(N))

# Simplified notation for sum calculation

sum_x = x[:]

The sum operation in JijModeling can also be performed on sets. In the combinatorial optimization problem.

mathematical models often appear that take the sum over some set $V$, such as $$ \sum_{v \in V} v x_v $$. This implementation is as follows.

`# Define variables to be used`

V = jm.Placeholder('V', dim=1)

num_V = V.shape[0].

x = jm.Binary('x', shape=(num_V))

v = jm.Element('v', V)

# Compute the sum over a set V

sum_v = jm.Sum(v, v*x[v])

By defining `v`

as an element of a set in `jm.Element`

and writing `jm.Sum`

with `v`

as its first argument, it is easy to write an implementation of summing over a set `V`

. All that remains is to prepare concrete values of that set `V`

as instances.

### Multiple Sum Operations

To compute multiple sums, one may write them in a formula as follows.

This can be implemented in JijModeling as follows.

`# Define variables to be used`

Q = jm.Placeholder('Q', dim=2)

I = Q.shape[0].

J = Q.shape[1]

x = jm.Binary('x', shape=(I, J))

i = jm.Element('i', (0, I))

j = jm.Element('j', (0, J))

# Calculate the sum over the two indices i, j

sum_ij = jm.Sum([i, j], Q[i, j]*x[i, j])

If there are multiple sums, multiple `jm.Sum`

can be omitted by making the subscripts a list `[subscript1, subscript2, ...]`

.
Of course, this can be written as $\sum_{i, j} = \sum_{i} \sum_{j}$, as can be seen from

`sum_ij = jm.Sum(i, jm.Sum(j, Q[i, j]*x[i, j]))`

### Sum Operation with condition

In addition to summing over the entire range of possible indexes, you may want to sum only if the indexes satisfy certain conditions, as in the following.

This can be implemented using JijModeling as follows

`# Define variables to be used`

I = jm.Placeholder('I')

x = jm.Binary('x', shape=(I))

i = jm.Element('i', (0, I))

U = jm.Placeholder('U')

# Calculate sum only for parts satisfying i<U

sum_i = jm.Sum((i, i<U), x[i])

In the part of `jm.Sum`

specifying the index, use a tuple to specify `(index, condition)`

like this.

The conditional expressions satisfied by the subscripts can be `<, <=, >=, >, ==, ! =`

can be used for the conditional expression that the subscript satisfies.

If you want to impose multiple conditions on a subscript, for example.

Suppose you want to implement a formula such as This can be implemented as follows

`# Define variables to be used`

d = jm.Placeholder('d', dim=1)

I = d.shape[0].

x = jm.Binary('x', shape=(I))

U = jm.Placeholder('U')

N = jm.Placeholder('N')

i = jm.Element('i', (0, I))

# Calculate sum only for the part satisfying i<U and i≠N

sum_i = jm.Sum((i, (i<U)&(i!=N)), d[i]*x[i])

When multiple conditions are imposed, the operator `&`

for logical AND and `|`

for logical OR can be used to chain the conditions. Of course, more than three conditions can be imposed.

If there is a condition on the subscripts in a sum operation on multiple subscripts, e.g..

A formula such as the following can be implemented as follows

`# Define variables to be used`

R = jm.Placeholder('R', dim=2)

I = R.shape[0]

J = R.shape[1]

x = jm.Binary('x', shape=(I, J))

i = jm.Element('i', (0, I))

j = jm.Element('j', (0, J))

N = jm.Placeholder('N')

L = jm.Placeholder('L')

# Calculate sum only for the part satisfying i>L and i!=N and j<i

sum_ij = jm.Sum([(i, (i>L)&(i!=N)), (j, j<i)], R[i, j]*x[i, j])

The part of `jm.Sum`

specifying the indexes should be replaced with `[[(index 1, condition of index 1), (index 2, condition of index 2), ...]] By making it look like `

[(index1, condition of index1), (index2, condition of index2), ...]', multiple sum operations can be written with conditions attached to each index.

[(j, j<i), (i, (i>L)&(i!=N))] will occur an error because $i$ is not yet defined in the $j<i$ part. This can be written in a formula like $\sum_{\substack{i>L \\ i!=N}} \left( \sum_{j<i} \cdots \right)$, but $\sum_{j<i} \left( \sum_{\substack{i>L \\ i!=N}} \cdots \right)$ corresponds to the fact that it cannot be written as Please note the order in which conditions are imposed on subscripts in multiple sums.

### Notes on Replacing Sums by Variables

For example, consider the following implementation of the formula.

Since $\sum_{j>i} x_{ij}$ appears twice, the following can be implemented by temporarily replacing it with an appropriate variable.

`# Define the variable to be used`

N = jm.Placeholder('N')

x = jm.Binary('x', shape=(N, N))

i = jm.Element('i', (0, N))

j = jm.Element('j', (0, N))

# Calculate sum only for j satisfying j>i

sum_j = jm.Sum((j, j>i), x[i, j])

# Compute sum for i

H = jm.Sum(i, sum_j*(sum_j-1))

Note that if `jm.Sum`

is replaced by a variable, it will not be written as `sum_j[i]`

in the following.

## Adding constraints

### Equality Constraints

Let us insert the following constraint condition into the problem.

By using an abbreviated description of the sum, we can easily implement the following.

`# Define variables to be used`

N = jm.Placeholder('N')

x = jm.Binary('x', shape=(N))

# Create a problem

problem = jm.Problem('constraint_example_01')

# Add constraint: Σ x = 1

problem += jm.Constraint('equality', x[:]==1)

Set a constraint with `jm.Constraint`

. The first argument is a string that is the name of that constraint, and the second argument is an equation with `==`

representing an equality constraint.

The constraint so created is then added to `problem`

using the `+=`

operator.

### Inequality Constraints

Using `jm.Constraint`

and `<=`

from earlier, you can also easily implement inequality constraints. For example, consider implementing the following equation.

The implementation using JijModeling is as follows.

`# Define variables to be used`

N = jm.Placeholder('N')

x = jm.Binary('x', shape=(N))

y = jm.Binary('y', shape=(N))

U = jm.Placeholder('U')

# Create problem

problem = jm.Problem('constraint_example_02')

# Add constraint: Σ x < Σ y + U

problem += jm.Constraint('inequality', x[:]-y[:]<=U-1)

We can only use `>=, <=`

as a constraint in JijModeling.

### Using forall to define constraints in a batch

There are many combinatorial optimization problems where there are $N$ constraints for subscripts as follows.

Such constraints can also be implemented cleanly in JijModeling.

`# Define variables to be used`

N = jm.Placeholder('N')

x = jm.Binary('x', shape=(N, N))

i = jm.Element('i', (0, N))

# Create problem

problem = jm.Problem('constraint_example_03')

# Add constraint: Σ_j x_{ij} = 1 (∀ i)

problem += jm.Constraint('forall', x[i, :]==1, forall=i)

By entering a subscript defined in `jm.Element`

as the `forall`

argument of `jm.Constraint`

, you can set constraints on all possible ranges for that subscript at once.

### Multiple forall

As the number of subscripts increases, constraint conditions may be imposed on multiple $\forall$ symbols as follows.

Similar to the way `jm.Sum`

sums over multiple subscripts, `[subscript 1, subscript 2, ...]]`

, this can be implemented by listing them as follows.

`# Define variables to be used`

N = jm.Placeholder('N')

x = jm.Binary('x', shape=(N, N, N))

i = jm.Element('i', (0, N))

k = jm.Element('k', (0, N))

# Create problem

problem = jm.Problem('constraint_example_04')

# Add constraint: Σ_j x_{ijk} = 1 (∀ i, k)

problem += jm.Constraint('multi_forall', x[i, :, k]==1, forall=[i, k])

### Attaching conditions to forall

Unlike before, we may wish to impose a constraint on all subscripts satisfying a certain condition as follows.

As with `jm.Sum`

, this can also be implemented by entering a tuple like `(subscript, condition)`

.

`# Define variables to be used`

A = jm.Placeholder('A', dim=1)

N = A.shape[0].

x = jm.Binary('x', shape=(N, N, N))

i = jm.Element('i', (0, N))

j = jm.Element('j', (0, N))

k = jm.Element('k', (0, N))

# Create problem

problem = jm.Problem('constraint_example_05')

# Add constraint: Σ_i x_{ijk} = A_k (∀ j, k | j > k)

problem += jm.Constraint('multi_forall_condition', x[:, j, k]==A[k], forall=[k, (j, j>k)])

Also, if you want to set conditions on `forall`

for multiple indexes, you can use `[[(index1, condition for index1), (index2, condition for index2), ...]]`

.
It is also possible to impose a condition on each index by doing. Multiple conditions on a single index can be written using the AND operator `&`

and the OR operator `|`

.

A statement like [(j, j>k), k] will occur an error because $k$ is not yet defined in the $j>k$ part.
Please be careful about the order of conditions in `forall`

for multiple subscripts.

### Variable Fixation by Obvious Constraints

When mathematically modeling a combinatorial optimization problem, it is sometimes possible to fix variables based on trivial constraints. In JijModeling, this can be easily implemented by using `jm.Constraint`

. For example, to fix $x_{0, 0}=0$ in a binary variable sequence $x_{i, j}$, use the following

`# Define the variables to be used`

N = jm.Placeholder('N')

x = jm.Binary('x', shape=(N, N))

# Create problem

problem = jm.Problem('constraint_example_06')

# Add an obvious constraint: x_{0, 0} = 1

problem += jm.Constraint('single_fix', x[0, 0]==1)

When such an obvious constraint is added, $x_{0, 0}$ is not passed to the solver as a decision variable, but is translated into an instance fixed at 0.

However, the obvious constraints must satisfy the following requirements.

- It is an equational constraint
- The left-hand side must be a decision variable only and must not contain other operations such as
`jm.Sum`

. - The right-hand side must not contain a decision variable

It is also possible to set them all at once using `forall`

. For example,

would be implemented as follows.

`# Define variables to be used`

N = jm.Placeholder('N')

x = jm.Binary('x', shape=(N, N))

j = jm.Element('j', (0, N))

# Create problem

problem = jm.Problem('constraint_example_07')

# Add an obvious constraint: x_{0, j} = 0 (∀ j)

problem += jm.Constraint('forall_fix', x[0, j]==1, forall=j)

## Check the implementation

If you are using Jupyter Notebook, you can check the mathematical model implemented in JijModeling while implementing. The following figure shows the actual implementation of equations (1)~(3) on Jupyter Notebook.

By using `.set_latex`

in the variable definition, you can set the formula display on Jupyter Notebook. If not set using `.set_latex`

, the variable `N`

defined in `d.shape[0]`

will be displayed as $d_mathrm{shape(0)}$ by default, as follows.

## Solving with JijZept

JijZept subscribers can solve optimization problems by sending the mathematical model implemented in JijModeling to the JijZept server. In doing so, the numerical data of the actual instance to be solved is also passed at the same time.

For example, the solution to the previous implementation is as follows

`import jijzept as jz`

# Prepare the instance

instance_data = {'d': [1, 10, -5, 0.45, -0.1]}

# Set up sampler

sampler = jz.JijSASampler(token='*** your API key ***', url='https://api.jijzept.com')

# Calculate using the sampler

response = sampler.sample_model(problem, instance_data, {'onehot': 5.0})

If a constraint is inserted into `problem`

using `jm.Constraint`

, the third argument to `sampler.sample_model`

is the undetermined multiplier of that constraint, of type `{'constraint name': undetermined multiplier number}`

in the dictionary. If there is no constraint, specify an empty dictionary.

JijZept uses the `.sample_model`

of the `Sampler`

class to submit the problem to be solved. See the JijZept documentation for more details.

## Solving with OpenJij in the local environment

You can also solve optimization problems implemented in JijModeling in your local environment by converting the JijModeling implementation to PyQUBO, which can then be used in OpenJij PyQUBO and other solvers.

The conversion to PyQUBO uses the information of a previously prepared instance.

`import openjij as oj`

# Convert to PyQUBO

pyq_obj = problem.to_pyqubo(instance_data, {})

# Convert for OpenJij using PyQUBO functionality

bqm = pyq_obj.compile().to_bqm(feed_dict={'onehot': 5.0})

# Prepare SA sampler

sampler = oj.SASampler()

# Solve using SA

response = sampler.sample(bqm)

## Decoding the computation, viewing the results

Both JijZept and OpenJij return the result in the form of a `response`

. However, this `response`

has a structure that is difficult to see when checking the obtained solution. Therefore, JijModeling decodes the `response`

to analyze and visualize whether the obtained results satisfy the constraint conditions.

The decoding is done as follows

`decoded = problem.decode(response, instance_data, {})`

Decoding is done by passing the result, `response`

, and the instance, `instance_data`

, to `.decode`

of the `Problem`

class.

Also, to check the state of the solution (decision variable sequence), the value of the objective function, and the constraint breaking from the `decode`

obtained in this way, do the following.

`print(decoded.solutions) # <- check solution state`

print(decoded.objectives) # <- check the value of the objective function

print(decoded.constraint_violations) # <- check for constraint violation

Running this, for example, yields the following output for an optimization problem using the mathematical model described in equations (1)~(3) above.

`[{'x': array([0., 0., 1., 0., 0.])}]`

[-5.]]

[{'onehot': 0.0}]

`decoded.solutions`

gives the state of the resulting computation. In this case, the decision variable sequence`x`

is used, and what kind of {0, 1} number sequence it became is converted to`numpy.ndarray`

and output.- In
`decoded.objects`

, the value of the objective function obtained from the result of the obtained calculation is obtained. In this case, the value -12.0 is obtained. - In
`decoded.constraint_violations`

, you can see the constraint violations from the obtained calculation results. Here only one constraint,`onehot`

, is shown. If multiple constraints are added, the degree to which each constraint is violated is displayed.

It is also possible to extract only those solutions for which all constraints are satisfied (feasible solutions). To do so, use the following

`print(decoded.feasibles().solution) # <- check the status of feasible solutions`

print(decoded.feasibles().objective) # <- check the value of the objective function of the feasible solution