Abstract: Naturally fractured reservoirs make important contributions to global oil and gas reserves and production. The modeling and simulation of naturally fractured reservoirs are different from conventional reservoirs as the existence of natural fractures. To address the development optimization problem of naturally fractured reservoirs, we propose an optimization workflow by coupling the optimization methods with the embedded discrete fracture model (EDFM). Firstly, the effective and superior performance of the workflow is verified based on the conceptual model. The stochastic simplex approximate gradient (StoSAG) algorithm, the ensemble optimization (EnOpt) algorithm, and the particle swarm optimization (PSO) algorithm are implemented for the production optimization of naturally fractured reservoirs based on the improved versions of the Egg model and the PUNQ-S3 model. The results of the two cases demonstrate the effectiveness of this optimization workflow by finding the optimal well controls which yield the maximum net present value (NPV). Compared to the initial well control guess, the final NPV obtained from the production optimization of fractured reservoirs based on all three optimization algorithms is significantly enhanced. Compared with the optimization results of the PSO algorithm, StoSAG and EnOpt have significant advantages in terms of final NPV and computational efficiency. The results also show that fractures have a significant impact on reservoir production. The economic efficiency of fractured reservoir development can be significantly improved by the optimization workflow.