Thanks for your explanations, Kris and Park.
Kris: I created 100,000 folders and each folder contained an individual case. The initial and boundary conditions are all different for the cases, so the Stochastic mode does not work in my case.
Park: I used LSF for the request. All scripts are attached. bsub_1.sh is the script submitted to LSF. The python script ctr_1.py takes the cases number and then call the script pflotran.sh to run all the cases. It worked well in a small test that requested 400 cores with one case running in each core. However, when the number of requested core increased to 1000, many cores just did nothing. I don't know if it is because the job submissions was delayed if the number of pending jobs reached the limit. I received the message "User <ese-chenkw>: Pending job threshold reached. Retrying in 60 seconds..." when I submit too many jobs at a time. Thanks.