The development of a new hard x-ray beamline I-TOMCAT equipped with a 1-meter-long short-period bulk high-temperature superconductor undulator (BHTSU) has been scheduled for the upgrade of the Swiss Light Source (SLS 2.0) at the Paul Scherrer Institute (PSI). The very hard x-ray source generated by the BHTSU will increase the brilliance at the beamline by over one order of magnitude in comparison to other state-of-the-art undulator technologies and allow experiments to be carried out with photon energies in excess of 60 keV. One of the key challenges for designing a 1-meter-long (100 periods) BHTSU is the large-scale simulation of the magnetization currents inside 200 staggered-array bulk superconductors. A feasible approach to simplify the electromagnetic model is to retain five periods from both ends of the 1-meter-long BHTSU, reducing the number of degrees of freedom (DOFs) to the scale of millions. In this paper, the theory of the recently-proposed 2D A-V formulation-based backward computation method is extended to calculate the critical state magnetization currents in the ten-period staggered-array BHTSU in 3D. The simulation results of the magnetization currents and the associated undulator field along the electron beam axis are compared with the well-known 3D H-formulation and the highly efficient 3D H-{\phi} formulation method, all methods showing excellent agreement with each other as well as with experimental results. The mixed H-{\phi} formulation avoids computing the eddy currents in the air subdomain and is significantly faster than the full H-formulation method, but is slower in comparison to the A-V formulation-based backward computation. Finally, the fastest and the most efficient A-V formulation in ANSYS 2020R1 Academic is adopted to optimize the integrals of the undulator field along the electron beam axis by optimizing the sizes of the end bulks.