Design of multiple branches splitting of equal mass flow rate in complex rheological flows like ice cream near melting point temperature can be a challenging task. Pulsations in flow rate due to air pumping process and small fluctuations in temperature affecting flow rheology can determine a consistent difference in internal pipe velocity distribution, resulting in a significant difference in the distribution of ice cream dosage. Computational sciences and engineering techniques have allowed a major change in the way products and equipment can be engineered, as a computational model simulating physical processes can be more easily obtained, rather than making prototypes and performing multiple experiments. Among such techniques, optimal shape design (OSD) represents an interesting approach. In OSD, the essential element with respect to classical numerical simulations in fixed geometrical configurations relays on the introduction a certain amount of geometrical degrees of freedom as a part of the unknowns. This implies that the geometry is not completely defined, but part of it is allowed to move dynamically in order to minimize or maximize an objective function. From a mathematical point of view, OSD is a branch of differentiable optimization and more precisely the application of optimal control for distributed systems. OSD is still today numerically difficult to implement, because it relies on a computer intensive activity and moreover because the concept of “optimal” is a compromise between shapes that are good with respect to several criteria. In this work, the applications of a multivariate constrained optimization algorithm is proposed in the case of a mechanical ice cream 1 to 5 splitting system, required to distribute in an evenly way from one freezer into five dosing valves. Results allowed to design a retro-fitting system on an existing production plant reducing the dosing error down to 3% on the average.