## 1 Introduction  
Topology optimization of multi-materials is an extension of single-material-based topology optimization and can be used in the structural optimization of components, modules, and also mechanical systems, which has become a hot topic in recent years [1; 2]. The challenge of topology optimization of multi-materials is exploring the appropriate mathematical formulation of the interpolation function and the optimization model, which enable each material in the design domain to be completely covered without overlapping, thus being easily manufactured [3; 4; 5; 6]. The density-based method is one of the widely used methods for topology optimization of multi-materials, which extends Solid Isotropic Material with Penalization (SIMP), Discrete Material Optimization (DMO), and other interpolation schemes for single materials [7; 8]. [PERSON] and [PERSON] [9] first extended the SIMP interpolation method to the multi-material problem and obtained microstructures with extreme thermal expansion. [PERSON] and [PERSON] [10] presented two multi-material interpolation formats that include an arbitrary number of material phases based on the SIMP and Rational Approximation of Material Properties (RAMP) methods [7]. In addition, [PERSON] and [PERSON] [11] proposed a discrete material optimization method for the topology optimization of composite structures, which was first adopted for the topology optimization of multi-materials by [PERSON] and [PERSON] [12]. [PERSON] [13] proposed a shape function parameterization (SFP) method to reduce the number of design variables in DMO. [PERSON] et al. [7] proposed a multi-material formulation using the ground structure method (GSM) considering material nonlinearity and further derived a ZPR update scheme for topology optimization of truss structures, which is effective for handling multiple constraints of volume/mass. Both the SIMP- and the DMO-based methods used multiple design variables to represent multi-materials, while [PERSON] and [PERSON] [14] proposed a peak function that used a linear combination of a normal distribution function for the interpolation function, which only needs one design variable for the topology optimization of multi-materials. [PERSON] and [PERSON] [15] proposed an ordered multi-material SIMP interpolation function that only used a single design variable, and the corresponding Matlab code was given. As most of the conventional methods have the problem of gray elements for the optimized structure, especially at the interface between different materials, [PERSON] [16] proposed a mapping-based interpolation function for topology optimization of multi-materials to address these issues, which can obtain precise 0-1 results for each material. Although many research results have been published for the topology optimization of multi-materials, it is still difficult for new students and researchers to follow the work. An example is the topology optimization of single materials, which was first introduced by [PERSON] and [PERSON] [17] and has been widely used since the 99-line Matlab code was published by [PERSON] [18]. There have been many educational papers on single-material topology optimization for different methods such as the SIMP method, the bi-directional evolutionary structural optimization method (BESO), the level set method, the moving morphable components (MMC) method, the homogenization method, and so on. For the SIMP method, [PERSON] et al. [19] presented a compact 88-line Matlab code based on the 99-line Matlab code, which improved the calculation speed through matrix operation and also improved filtering strategies. [PERSON] and [PERSON] [20] presented an efficient and compact 169-line Matlab code to solve three-dimensional topology optimization problems. Later, [PERSON] [21] presented new topology optimization Matlab codes for 2D and 3D minimum compliance problems, consisting of 99 and 125 lines, respectively.