How to model and solve such a 0-1 programming problem

Question

My problem is described in this picture(It's like a Pyramid structure):

The objective function is below: $$\min\sum_{k=1}^\ell\sum_{i=0}^{2^k-1}\sum_{j=0}^{2^k-1}\left(A_{i,j}^k-A_{i,j+1}^k\cdot\frac{D_{i,j}^k}{D_{i,j+1}^k}\right)^2$$

In this problem, $D$ is known, $A$ is the object that I want to get. It is a layered structure, each block in the upper layer is divided into four sub-blocks in the layer below. And the value of the upper layer node is equal to the sum of the four child nodes of the lower layer. In above example, I used only 2 layers.

What I want to do is simulate the distribution of $D$ with $A$, so in the objective function is the ratio of two adjacent squares in each row in $A$ compared to the value in $D$. I do this comparison on each layer and sum them. Then it is all of my objective function. But in the finest layer, the value in $A$ is $0$ or $1$ like the above picture, so this is also a 0-1 programming problem. I have tried to solve it using integer quadratic programming in python library CVXPY but adding the constraints of A<=1. However, it seems the speed is slow.

So I want to model this problem as a standard 0-1 programming problem, and then use some solutions of 0-1 programming to solve it, but I don't know how to model this problem as a standard 0-1 programming problem, and what methods I can use to solve it?

In other way, if I use quadratic programming instead of Integer quadratic programming (from $A$ is $0/1$ to $A$ is a value between $0$ and $1$). In this way whether I can solve it using the method of derivation, because this is a convex optimization problem, which can guarantee the global optimal solution. There are two unknown variables in each item, that is, the two items with $A$ in the formula. Partial derivatives are obtained for them, and the restriction of A<=1 is added, then solve using (such as) gradient descent method. Is this mathematically feasible, because I don't know much about optimization, and if it is possible, how should I do it? If not possible, what other methods can I use?

Answer 1

Say the finest layer, $k$ is 1, then declare $n \times m$ variables $a_{n,m}^k \in\ \{0,1\}$ as binary for the finest layer. For other layers, $k=2...K$
$A_{i,j}^{k+1} = a_{i,j}^k+a_{i+1,j+1}^k+a_{i+1,j}^k+a_{i,j+1}^k$.

Except the finest layer $k=1$, rest of $a$'s are real numbers.